Computer Math II

UTG

CPE 332

Computer Engineering

Mathematics II

Chapter 1 Vector

Web Site

http://cpe.rsu.ac.th/ut

– Download Material, Course Notes

– Download Slides

– Download HW Solutions

– Grading

– Announcements

– Resources

Today Topics

• Period 1

– Course Outlines

– Course Web Site

– Part I Chapter 1 Vector (Review)

– Breaks

– Part I Chapter 1 Vector (Review)

• Assignment:

– ยังไม่มีการบ ้าน

Download MATLAB Tutorial 1-5 และลองทํา Exercise ดู

• Next Week ต่อ Vector และ Chapter 2 เรื่อง Matrix

CPE 332 T1-57 Wk1

Definition of Vector

Definition of Vector

Notes

• เนื่องจาก Vector มีทั้งขนาดและทิศทาง เรา

สามารถเขียน Vector เป็นสองส่วน

– ส่วนขนาดแทนที่ด ้วย Scalar

– ส่วนทิศทาง จะแทนที่ด ้วย Unit Vector ที่มีทิศทาง

เดียวกับ Vector เดิม

ˆ

ˆ

F = F = F f = F f = F i + F j + F k = F (cosαi + cos βj + cos k γ )

1

2

3

• การกําหนดทิศทาง อาจจะกําหนดเป็น Component ในแกน

Coordinate (x,y,z); อาจจะกําหนดเป็นมุมที่กระทํากับแกน

Coordinate

• อาจจะกําหนดเป็น Ratio ที่กระทํากับแกนก็ได ้

• จะกล่าวต่อไปภายหลัง

– เราจะเน ้นที่สองอันแรก คือกําหนดเป็น Component i,j,k ในแกน

x,y,z

– หรือกําหนดในรูป Cosine ของมุม

– ทั้งสองอันนี้จะเกี่ยวข ้องกับ Unit Vector

Vector Operations

• เนื่องจาก Vector ประกอบด ้วยทั้งขนาดและ

ทิศทาง

– พีชคณิต เช่น บวก ลบ คูณ หาร จะไม่เหมือนกับ Scalar

เนื่องจากต ้องนําทิศทางมาประกอบการคํานวณด ้วย

– การ บวก-ลบ ของ Vector จะได ้ Vector ใหม่ที่ขนาด

และทิศทางต่างจากเดิม

– การคูณ เราจะไม่ใช ้คําว่า ‘Multiplication’ แต่จะใช ้คําว่า

‘Product’ แบ่งเป็นสองประเภท

• Scalar Product (Dot Product; ●) จะได ้ Scalar

• Vector Product (Cross Product; X) จะได ้ Vector ที่ตั้งฉาก

กับ Vector เดิมทั้งสอง

Addition and Substraction

การประยุกต์ใช ้ใน Plane

Geometry

Component Vector

Component Vector in

Cartesian Coordinate

Component Vector in

Cartesian Coordinate

Component Vector in

Cartesian Coordinate

Position Vector

• จุดใน Space สามารถแสดงได ้โดยใช ้ Vector

เริ่มจาก Origin

– อาจเรียก Location Vector หรือ Radius Vector

– จุด P แสดงได ้โดยใช ้ Vector OP

– และสามารถแสดงได ้โดยใช ้ Component Vector

Position Vector และ

Addition-Subtraction using

Component Vector

สรุป

• การเขียน Vector ในลักษณะ Component

จะสามารถบวกและลบกันได ้ง่าย โดยการ

บวกลบแต่ละ Component บนแกนเดียวกัน

– Vector Product สามารถคํานวณได ้เช่นกัน

• จุดใน Space สามารถแทนด ้วย Vector เริ่ม

จากจุด Origin เรียก Position Vector

– Vector ที่เกิดจากสองจุดใน Space สามารถ

คํานวณได ้จาก Position Vector นี้

r

Any vectors in Cartesian

Coordinates

• Given 2 Points, P(x1,y1,z1) and

Q(x2,y2,z2)

– We have OP+PQ=OQ

– Then PQ = OQ – OP

• PQ = x2i+y2j+z2k – x1i+y1j+z1k

• PQ =(x2-x1)i+(y2-y1)j+(z2-z1)k

Z

Q(x2,y2,z2)

O

Y

P(x1,y1,z1)

X

Any vectors in Cartesian

Coordinates

• Given 2 Points, P(x1,y1,z1) and

Q(x2,y2,z2)

– PQ =(x2-x1)i+(y2-y1)j+(z2-z1)k

– Also magnitude or length of vector is the

distance between those 2 points (Euclidian

Distance)

• PQ = √(x2-x1)2+(y2-y1)2+(z2-z1)2

Z

Q(x2,y2,z2)

O

Y

P(x1,y1,z1)

X

Direction Cosine/Ratio

• Vector สามารถเขียนเป็นสองส่วนประกอบ

=

A

A aˆ

– ขนาด สามารถหาได ้ง่าย กรณี Position Vector

– ทิศทาง คือ Unit Vector ที่มีทิศทางเดียวกันกับ

Vector นั้น

• ทิศทาง สามารถแตกเป็น Component Vector บนแต่ละ

แกนได ้ด ้วย

• ทิศทางสามารถกําหนดด ้วยมุมที่ทํากับแต่ละแกนได ้ด ้วย

• ทั้งสองแบบนี้ สัมพันธ์กันทางตรีโกณมิติ โดยการกําหนด

ด ้วยค่า Cosine ของมุม เรียก Direction Cosine

Direction Cosine

• Position vector OP

– Magnitude equal to OP = √x2+y2+z2

– Direction: cosαi+cosβj+cosγk

Called Direction Cosine

We have

F3

cosα =F1/OP

cosβ =F2/OP

F2

cosγ =F3/OP

F1

Direction Cosine and

Direction Ratio

Direction Cosine and

Direction Ratio

Example

• Given points P1(2,-4,5) and P2(1,3,-2),

find the vector P1P2 and its magnitude

and direction

– OP1 = 2i-4j+5k and OP2 = i+3j-2k

– P1P2=OP2-OP1=-i+7j-7k

– P1P2 = √1+49+49=√99

– Cos α = -1/√99 then α = 95.8 degree

– Cos β = 7/√99 then β = 45.3 degree

– Cos γ = -7/√99 then γ = 134.7 degree

Direction Cosine and

Direction Ratio

Next Week

• Vector Product

– Scalar Product(Dot)

– Vector Product(Cross)

• Chapter II: MATRICES

• HW II

CPE 332

Computer Engineering

Mathematics II

Week 2

Chapter 1 Vector (cont.)

Chapter 2 Matrix

Today Topics

• Chapter 1 Cont.

• Break

• Chapter 2: Matrix

• Download Homework 1: Chapter 1

– Due Next Week

Component Vector

Component Vector in

Cartesian Coordinate

Component Vector in

Cartesian Coordinate

Position Vector

• จุดใน Space สามารถแสดงได ้โดยใช ้ Vector

เริ่มจาก Origin

– อาจเรียก Location Vector หรือ Radius Vector

– จุด P แสดงได ้โดยใช ้ Vector OP

– และสามารถแสดงได ้โดยใช ้ Component Vector

Position Vector และ

Addition-Subtraction using

Component Vector

Any vectors in Cartesian

Coordinates

• Given 2 Points, P(x1,y1,z1) and

Q(x2,y2,z2)

– PQ =(x2-x1)i+(y2-y1)j+(z2-z1)k

– Also magnitude or length of vector is the

distance between those 2 points (Euclidian

Distance)

• PQ = √(x2-x1)2+(y2-y1)2+(z2-z1)2

Z

Q(x2,y2,z2)

O

Y

P(x1,y1,z1)

X

Direction Cosine/Ratio

• Vector สามารถเขียนเป็นสองส่วนประกอบ

=

A

A aˆ

– ขนาด สามารถหาได ้ง่าย กรณี Position Vector

– ทิศทาง คือ Unit Vector ที่มีทิศทางเดียวกันกับ

Vector นั้น

• ทิศทาง สามารถแตกเป็น Component Vector บนแต่ละ

แกนได ้ด ้วย

• ทิศทางสามารถกําหนดด ้วยมุมที่ทํากับแต่ละแกนได ้ด ้วย

• ทั้งสองแบบนี้ สัมพันธ์กันทางตรีโกณมิติ โดยการกําหนด

ด ้วยค่า Cosine ของมุม เรียก Direction Cosine

Direction Cosine

• Position vector OP

– Magnitude equal to OP = √x2+y2+z2

– Direction: cosαi+cosβj+cosγk

Called Direction Cosine

We have

F3

cosα =F1/OP

cosβ =F2/OP

F2

cosγ =F3/OP

F1

Products of Vectors

• Vector Product

– Scalar Product(DOT)

– Vector Product(Cross)

Scalar Product(DOT)

Scalar Product (DOT)

Scalar(Dot) Product

A

θ

n

A●n=Acosθ

Scalar(Dot) Product

A●(B+C)=A●B+A●C

– Let A = a1i+a2j+a3k, B = b1i+b2j+b3k

• We have A●B = a1b1+a2b2+a3b3

– Also

– Given S=ai+bj, the equation of line

perpendicular to this vector is in the form

ax+by=c

Line ax+by=c

S=ai+bj

DOT Product

Example

• Find the angle between the vector

A=i-j-k and B = 2i+j+2k

• We calculate A●B = 1.2-1.1-1.2=-1

• Also A = √(1+1+1)=√3

• Also B = √(4+1+4)=3

• Then Cos θ = -1/3√3

– θ = 101.1 degrees

Vector Product (Cross)

Cross Product

3 Vector Products

Examples

• Let A=2i+3j-k, B=i+j+2k

– A●B = 2+3-2 = 3

– A×B = (6+1)i-(4+1)j+(2-3)k=7i-5j-k

– A×B is orthogonal to both A and B

• Test : A●(A×B) = (2i+3j-k)●(7i-5j-k) = 14-

15+1=0

• Test : B●(A×B) = (i+j+2k)●(7i-5j-k) = 7-5-

2=0

Plane Equation in 3D

• ใน 2D สมการเส ้นตรงจะมี general Form

– Ax+By=C

• ใน 3D สมการของ Plane จะมี General

Form

– Ax+By+Cz=D

– D เป็นค่าคงที่ ทุกๆสมการในรูปเดียวกัน แต่ค่า

D ต่างกัน จะเป็นระนาบที่ขนานกัน

• 3x-2y+5z = 3 จะขนานกับ 3x-2y+5z = 6

Example 1

• กําหนดสมการของ Plane 2x+3y+2z=5 จง

หา unit vector ที่ตั้งฉากกับ Plane นี้

– กําหนด 3 จุด คือ A, B, C ดังนี้

– A: x=0,y=0,ดังนั้น z=5/2  A(0,0,2.5)

– B: x=1,y=0, ดังนั้น z=(5-2)/2  B(1,0,1.5)

– C: x=0,y=1, ดังนั้น z=(5-3)/2  C(0,1,1)

– Vector AB x AC จะได ้ Vector ที่ตั้งฉากกับ

Plane

AB = i k,

AC = j − 1. k

5

i

j

k

AB × AC = 1

0

−1 = i +1.5j + k

0

1

−1.5

1

uˆ =

[i +1 5. j+ k]

.

4 25

• สังเกตุว่าทุกๆ Vector ที่เป็น multiple ของ

2i+3j+2k จะตั้งฉากกับ Plane 2x+3y+2z=k

เสมอ โดยที่ k เป็นค่าคงที่ใดๆ

Example 2

• จงหาสมการของ Plane ที่ตั้งฉากกับ Vector

3i-2j-k และกําหนดให ้จุด (1,1,2) อยู่บน

Plane นั้น

– จากตัวอย่างก่อน เราได ้สมการของ Plane เป็น

3x-2y-z= k

– เราหาค่า k โดยแทนค่าจุด (1,1,2) ลงใน

สมการดังนี้ 3(1)-2(1)-(2)=-1=k

– ดังนั้นสมการที่ต ้องการจะเป็น 3x-2y-z+1=0

Definition of Matrix

Row Matrix, Column Matrix

Basic Operations

Matrix Multiplication

Square Matrix

Matrix Transpose

Types of Matrix

Types of Matrix

Types of Matrix

Types of Matrix

Types of Matrix

Matrix Inverse

Orthogonal/Unitary Matrix

Orthogonal Vector

Orthogonal Vector

Examples

1

3

5

A = 3

2

6

a

is

s

quare ma

trix of

orde

3

r ×

3, A

al

is

so a

s

ymmetric ma

trix

5

6

4

 0

− 3

5 

B =

3

0

− 6

a

is

s

quare ma

trix of

orde

3

r ×

3, B

al

is

so a

sk

ew - symmetric ma

trix

− 5 6

0 

1

0

0

C = 0

2

0

a

is

di

agonal ma

trix, Di

ag( )

C = {1,

3

2,

t

}, race( )

C = 1 + 2 + 3 = 6

0 0 3

C

al

is

so a

s

ymmetric ma

trix.

Examples

2 3 3

D = 0

−1 4

an

is

uppe

di

r agonal ma

trix

0 0 1

2 0 0

E = 1

0

0

a

is

l

ower di

agonal ma

trix

3 0 1

1 0 0

1 0

I =

,

I =

2

3

0 1 0

0 1

0 0 1

Examples

F

a

or

ma

ny

trix A

ha

we

ve AAT

a

is

s

ymmetric ma

trix

, ATA

al

is

so a

s

ymmetric ma

trix.

T

T

T

T T

T T

T

T

proof : From (

AB) = B A

ha

we

ve (

AA ) = (A ) A = AA

0 0 0

F = 0

0

0

a

is

zer

o ma

trix of

orde

3

r × 3

0 0 0

0

G = [

 

0

0

0] a

is

zer

o row

vect

an

or

d H = 0

a

is

zer

o c

olumn ve t

c or

 

0

 

Examples

 3 

 

K = 3i − 2j + 4k can

w

be

ritten i

ma

n

trix

as

form

- 2

 4 

 

 2 

 

L = 2i + j k can

w

be

ritten i

ma

n

trix

as

form

1 

- 

 

1

 2 

T

 

K L = K L = [3 − 2 4] 1 = 3 ⋅ 2 + ( 2

− ) ⋅1+ (

4 − )

1 = 0

 

−1

 

in this

v

case

ector K an

v

d ector L ar

e ort

hogona (pe

l

rpendi ul

c ar.)

Examples

1 0 − 

1

 1

1

1 

 2

4

4

1

1

3 

M = 2

1

1 ,

N =

−

,

2

4

4 

0 1 − 

1

− 1 1

1 

 2 4

4 

1 0 − 

1  1

1

1 

1 0 0

  2

4

4

1

1

3 

ha

we

ve MN = 2

1

1

=

 

0

1

0

2

4

4 

0 1 − 

1 − 1

1

1 

0 0 1

  2 4

4 

can

We

say

that M =

1

-

N an

d N =

1

-

M si

n

ce MN = .

I

If ma

trix di

is

agonal t

, he i

nverse

is

just the ma

trix w h

it

reciprocal of

di

agonal el

ements an

d

al

is

so di

agonal.

1 0 0

1 0 0

−1

1

F

ex

or ampl

e, O = 0

2

0 then O

=

0

0

2

0 0 3

1

0

0

3 

Examples

 cosθ

sinθ 

P = 

− sinθ cosθ

an

is

ort

hogonal ma

trix si

nce

cosθ − sinθ 

P 1

-

= PT = 

sinθ

cosθ

can

we

,

t

see hat

 cosθ

sinθ cosθ

− sinθ 

PPT = 



− sinθ cosθ sinθ

cosθ 

cos2 θ + sin2 θ

− cosθ sinθ + sinθ cosθ 

= 

− sinθ cosθ + cosθ sinθ

c cos2 θ + sin2 θ

1

0

=

a

for

a

ny ngleθ

0 1

Determinant

Determinant of Matrix

• Given square matrix, determinant of a

matrix A written |A| is defined by

recursive equation as

n

i+

det(A) = A = ∑(− j

)

1

a Minor( a )

i, j

i, j

j =1

n

=

j +

(−

l

)

1

a

Minor( a )

j , l

j , l

j =1

– Starting from [a1x1] = a, and det(a)=a

Sign Matrix

• Given matrix A of order nxn

– Sign matrix of A is the matrix in the form

B=[bij]=[(-1)(i+j)]

a a a n

n

11

12

1 

 (− )

1 2

(− )

1 3

(−

+

)

1

1 

3

4

n + 2 

a

a

a n

21

22

2 

 (− )

1

(− )

1

(− )

1

A =  

 s

, ign A =

 

 

a

a

a

n +1

n + 2

2 n

n

n

nn

1

2

(− )

1

(− )

1

(− )

1

 1

−1 1 

−1 1

−1 

=  1 −1 1 

 

 

 

1 

Minor aij

• Is the determinant of matrix A after

taken out row ith and jth column.

a

a

a

a

a

 1

2

− 3

4

5 

11

12

13

14

15

 

a

a

a

a

a

− 2

3

4

− 5

6

21

22

23

24

25

 

A =  a

a

a

a

a  =  3

− 4

5

6

− 7

31

32

33

34

35

 

a

a

a

a

a

4

5

− 6

7

8

 41

42

43

44

45 

a

a

a

a

a

 − 5

6

7

8

− 9

51

52

53

54

55

a

a

a

a

1

2

4

5

11

12

14

15

a

a

a

a

3

− 4 6 − 7

Minor( a )

31

32

34

35

=

=

23

a

a

a

a

4

5

7

8

41

42

44

45

a

a

a

a

− 5

6

8

− 9

51

52

54

55

Cofactor aij

• Cofactor aij is the minor aij with the sign

according to sign pattern

• Matrix of cofactor of A is the matrix B

which each element bij is the cofactor

aij

 5 6

4

6

4

5 

1 2 3

 8 9

7

9

7

8 

 2 3

1

3

1

2 

.

Ex A = 4

5

6 , Co(A) =

−

7 8 9

 8 9

7

9

7

8 

 2 3

1

3

1

2 

 5 6

4

6

4

5 

Determinant of 2x2 and 3x3

Calculation of Determinant using

Recursive Expansion

n

det(A) = A = ∑

i+

(−

j

)

1

a Minor( a )

i, j

i, j

j =1

n

=

j +

(−

l

)

1

a

Minor( a )

j , l

j , l

j =1

a

b

a = a,

= ad bc

c

d

a

b

c

e

f

d

f

d

e

d

e

f = a

b

+ c

j

k

i

k

i

j

i

j

k

= a( ek jf ) − b( dk if ) + c( dj ie)

= aek-afi bdk + bfi + cdj cei, first row expand d

f

a

c

a

c

= b

+ e

j

i

k

i

k

d

f

= b

− ( dk if ) + e( ak ic) − j( af dc)

= bdk

+ bfi + aek cei afj + cdj, second column expand

Calculation of Determinant using

Recursive Expansion

n

det(A) = A = ∑

i+

(−

j

)

1

a Minor( a )

i, j

i, j

j =1

Complexity = O(n!)

n

=

j +

(−

l

)

1

a

Minor( a )

j , l

j , l

j =1

a

a

a

a

11

12

13

14

a

a

a

a

21

22

23

24

a

a

a

a

31

32

33

34

a

a

a

a

41

42

43

44

a

a

a

a

a

a

a

a

a

a

a

a

22

23

24

21

23

24

21

22

24

21

22

23

= a a

a

a

a a

a

a

+ a a

a

a

a a

a

a

11

32

33

34

12

31

33

34

13

31

32

34

14

31

32

33

a

a

a

a

a

a

a

a

a

a

a

a

42

43

44

41

43

44

41

42

44

41

42

43

ความหมายของ Determinant

• เป็นค่า Scalar ที่สําคัญที่สุดที่บ่งบอก

คุณสมบัติของ Square Matrix

คุณสมบัติของ Determinant

คุณสมบัติของ Determinant

Calculation of Determinant

(Algorithm)

n

det(A) = A = ∑

i+

(−

j

)

1

a Minor( a )

i, j

i, j

Complexity = O(n!)

j =1

n

=

j +

(−

l

)

1

a

Minor( a )

j , l

j , l

j =1

a

a

a

a

11

12

13

14

ถ ้าเราบวกลบ Column เพ่ือให ้ Element ในแถว(หรือคอลัมน์)

a

a

a

a

21

22

23

24

ที่ต ้องการขยายเป็นศูนย์หมดยกเว ้น Element เดียว

a

a

a

a

เราจะลงเอยด ้วยการคํานวณหา Determinant ของ Matrix

31

32

33

34

ที่มีขนาดลดลงหนึ่ง

a

a

a

a

41

42

43

44

a

a

a

a

a

a

a

a

a

a

a

a

22

23

24

21

23

24

21

22

24

21

22

23

= a a

a

a

a a

a

a

+ a a

a

a

a a

a

a

11

32

33

34

12

31

33

34

13

31

32

34

14

31

32

33

a

a

a

a

a

a

a

a

a

a

a

a

42

43

44

41

43

44

41

42

44

41

42

43

การบวกลบดังกล่าวต ้องมีหลักการ มิฉะนั้นผลลัพธ์จะไม่ถูกต ้อง

เราจะใช ้คุณสมบัติข ้อ 9 และ 10 ของ Determinant เพื่อกระทําดังกล่าว

Algorithm การหา Determinant ที่

มีประสิทธิภาพ

• ใช ้คุณสมบัติข ้อ 10 ร่วมกับข ้อ 9 เพื่อสร ้างเป็น Algorithm

– 1. มองหา Element ใน Matrix ที่มีค่าเท่ากับ 1 ถ ้าหาไม่ได ้ เลือก

Element ใดก็ได ้ จากนั้นหารทั้งแถว หรือหารทั้ง Column ด ้วยค่า

ของ Element นั้นเพื่อทําให ้ค่าเป็น 1 ตัวเลขที่มาหารนั้นจะต ้อง

กลับนํามาคูณกับคําตอบที่ได ้ เป็นค่า Determinant ที่ต ้องการ

(คุณสมบัติข ้อ 9)

– 2. พิจารณาว่าจะ Expand แบบแถวหรือ Column ผ่าน Element ที่

เลือก จากนั้นกําจัด Element อื่นในแนวที่ Expand เป็นศูนย์ให ้

หมด(คุณสมบัติข ้อ 10) สมมุติเราเลือก Element a(x,y)

• ถ ้าจะ Expand แบบแถว ให ้บวกลบ Column อื่นกับ Column ที่

ผ่าน Element ที่เลือก เพื่อให ้ Element ในแถวที่จะ Expand

เป็นศูนย์ทั้งหมด ยกเว ้น Element ที่เลือก

– Col(j) ใหม่ = Col(j) เก่า – a(x,j)*Col(y); j = 1,2,..,n ยกเว ้น y

• ถ ้าจะ Expand แบบ Column ให ้บวกลบแถวอื่นกับแถวที่ผ่าน

Element ที่เลือก เพื่อให ้ Element ใน Column ที่จะ Expand

เป็นศูนย์ทั้งหมด ยกเว ้น Element ที่เลือก

– Row(i) ใหม่ = Row(i) เก่า – a(i,y)*Row(x); i = 1,2,..,n ยกเว ้น x

Algorithm การหา Determinant ที่

มีประสิทธิภาพ(ต่อ)

– 3.ทําการ Expand ตามสูตร เราจะลงเอยด ้วย

การหา Determinant ของ Matrix ที่มีขนาด

ลดลงหนึ่งเพียงครั้งเดียว

– 4. วิธีนี้สามารถทําเป็น Recursive เพื่อลดการ

หา Determinant ของ Matrix ขนาดใหญ่

เหลือแค่การหา Determinant ของ Matrix 2x2

หรือ 3x3

การหา Determinant

การหา Determinant

Inverse of Matrix

Inverse of Matrix

Inverse of Matrix

Inverse of Matrix

Inverse of Matrix

Matrix Inverse

• การคํานวณตามสมการที่กล่าวมา จะใช ้การ

คํานวณมากเกินไป

• หลัง Midterm เราจะมาดู Algorithm ที่มี

ประสิทธิภาพ เพื่อใช ้ในการหา Matrix

Inverse

– Gauss-Jordan Method

Matrix Norms

Norms

เรียก Spectral Radius ของ X

Homework 1:Vector

• Download HW 1 Question Sheet จาก

Website

– พิมพ์ Question Sheet ลงบนกระดาษ A4

– ทําการบ ้าน โดยแสดงวิธีทําและคําตอบในช่อง

ที่กําหนด ด ้วยลายมือนักศึกษา ห ้ามพิมพ์ ต ้อง

ใช ้กระดาษที่พิมพ์นี้ทําการบ ้าน

– เขียนชื่อที่หัวกระดาษ ส่งอาทิตย์หน ้าต ้นชั่วโมง

– ไม่รับงานที่ไม่เป็นไปตามที่กําหนด

CPE 332

Computer Engineering

Mathematics II

Week 3: Ch.2 Matrices

Continue(Linear Equations)

Ch.3 Eigenvector

Today Topics

• Part I Chapt. 2 Matrices

• Break

• Chapter 2: Linear Equations

• Homework 1: Due

• Homework 2 ส่งสัปดาห์หน ้า

Chapter 2: Cont

• Last Week

– Matrix Concept and Notations

– Types of Matrix

– Concept of Determinant and Inverse

• This Week

– More on Determinant and Matrix

– System of Linear Equations: Homogeneous

and Non-homogeneous

– Eigenvalue/Eigenvector Concept

Determinant

Calculation of Determinant

n

det(A) = A = ∑

i+

(−

j

)

1

a Minor( a )

i, j

i, j

j =1

n

=

j +

(−

l

)

1

a

Minor( a )

j , l

j , l

j =1

a

b

a = a,

= ad bc

c

d

a

b

c

e

f

d

f

d

e

d

e

f = a

b

+ c

j

k

i

k

i

j

i

j

k

= a( ek jf ) − b( dk if ) + c( dj ie)

= aek-afi bdk + bfi + cdj cei, first row expand d

f

a

c

a

c

= b

+ e

j

i

k

i

k

d

f

= b

− ( dk if ) + e( ak ic) − j( af dc)

= bdk

+ bfi + aek cei afj + cdj, second column expand

คุณสมบัติของ Determinant

คุณสมบัติของ Determinant

Calculation of Determinant

n

det(A) = A = ∑

i+

(−

j

)

1

a Minor( a )

i, j

i, j

Complexity = O(n!)

j =1

n

=

j +

(−

l

)

1

a

Minor( a )

j , l

j , l

j =1

a

a

a

a

11

12

13

14

ถ ้าเราบวกลบ Column เพ่ือให ้ Element ในแถว(หรือคอลัมน์)

a

a

a

a

21

22

23

24

ที่ต ้องการขยายเป็นศูนย์หมดยกเว ้น Element เดียว

a

a

a

a

เราจะลงเอยด ้วยการคํานวณหา Determinant ของ Matrix

31

32

33

34

ที่มีขนาดลดลงหนึ่ง

a

a

a

a

41

42

43

44

a

a

a

a

a

a

a

a

a

a

a

a

22

23

24

21

23

24

21

22

24

21

22

23

= a a

a

a

a a

a

a

+ a a

a

a

a a

a

a

11

32

33

34

12

31

33

34

13

31

32

34

14

31

32

33

a

a

a

a

a

a

a

a

a

a

a

a

42

43

44

41

43

44

41

42

44

41

42

43

การบวกลบดังกล่าวต ้องมีหลักการ มิฉะนั้นผลลัพธ์จะไม่ถูกต ้อง

เราจะใช ้คุณสมบัติข ้อ 9 และ 10 ของ Determinant เพื่อกระทําดังกล่าว

Algorithm การหา Determinant ที่

มีประสิทธิภาพ

• ใช ้คุณสมบัติข ้อ 10 ร่วมกับข ้อ 9 เพื่อสร ้างเป็น Algorithm

– 1. มองหา Element ใน Matrix ที่มีค่าเท่ากับ 1 ถ ้าหาไม่ได ้ เลือก

Element ใดก็ได ้ จากนั้นหารทั้งแถว หรือหารทั้ง Column ด ้วยค่า

ของ Element นั้นเพื่อทําให ้ค่าเป็น 1 ตัวเลขที่มาหารนั้นจะต ้อง

กลับนํามาคูณกับคําตอบที่ได ้ เป็นค่า Determinant ที่ต ้องการ

(คุณสมบัติข ้อ 9)

– 2. พิจารณาว่าจะ Expand แบบแถวหรือ Column ผ่าน Element ที่

เลือก จากนั้นกําจัด Element อื่นในแนวที่ Expand เป็นศูนย์ให ้

หมด(คุณสมบัติข ้อ 10) สมมุติเราเลือก Element a(x,y)

• ถ ้าจะ Expand แบบแถว ให ้บวกลบ Column อื่นกับ Column ที่

ผ่าน Element ที่เลือก เพื่อให ้ Element ในแถวที่จะ Expand

เป็นศูนย์ทั้งหมด ยกเว ้น Element ที่เลือก

– Col(j) ใหม่ = Col(j) เก่า – a(x,j)*Col(y); j = 1,2,..,n ยกเว ้น y

• ถ ้าจะ Expand แบบ Column ให ้บวกลบแถวอื่นกับแถวที่ผ่าน

Element ที่เลือก เพื่อให ้ Element ใน Column ที่จะ Expand

เป็นศูนย์ทั้งหมด ยกเว ้น Element ที่เลือก

– Row(i) ใหม่ = Row(i) เก่า – a(i,y)*Row(x); i = 1,2,..,n ยกเว ้น x

Algorithm การหา Determinant ที่

มีประสิทธิภาพ(ต่อ)

– 3.ทําการ Expand ตามสูตร เราจะลงเอยด ้วย

การหา Determinant ของ Matrix ที่มีขนาด

ลดลงหนึ่งเพียงครั้งเดียว

– 4. วิธีนี้สามารถทําเป็น Recursive เพื่อลดการ

หา Determinant ของ Matrix ขนาดใหญ่

เหลือแค่การหา Determinant ของ Matrix 2x2

หรือ 3x3

การหา Determinant

Column 3 Expansion

Or Row 2 Expansion

การหา Determinant

Inverse of Matrix

Inverse of Matrix

Inverse of Matrix

Inverse of Matrix

Inverse of Matrix

สาธิตทดลอง Row Expansion

Matrix Norms

Norms

เรียก Spectral Radius ของ X

Break

System of Linear Equations

System of Linear Equations

Reduced Matrix

Reduced Matrix

Reduced Matrix

Solutions of Homogeneous

Solutions of Homogeneous

Solutions of Homogeneous

Non-homogeneous Systems

Non-homogeneous Systems

Non-homogeneous Systems

Non-homogeneous Systems

Non-homogeneous Systems

Non-homogeneous Systems

Chapter 3 Intro

MATLAB TUTORIAL

End of Week 3

• Download HW 2 Due Next Week

– Chapter 1: Vector Operations

– Chapter 2: Linear Equations

• Next week: WK 4

– Eigenvalue/Eigenvector Continue

– MATLAB

– HW 3

CPE 332

Computer Engineering

Mathematics II

Week 4: Ch.3 Eigenvector and

Diagonalization

Ch4. Probability Review

Today Topics

• Part I Chapt. 3 Eigenvalue/Eigenvector

• Break

• Chapter 4: Review Probability

• Homework 2: Due

• Homework 3 ส่งสัปดาห์หน ้า ต ้นชั่วโมง

Determinant

Chapter 3 Intro

Diagonalization

Orthogonal and Symmetric

• จาก Theorem ที่กล่าว จุดที่สําคัญคือ

Theorem ที่ 8 สุดท ้าย กล่าวคือ

– ถ ้า Matrix เป็น Real Symmetric เราจะได ้

Eigenvalue เป็นค่าจริง และจะมี Orthogonal

Matrix ที่จะมา Diagonalize ได ้

– เมื่อ Matrix เป็น Orthogonal มันจะหา Inverse

ได ้ง่าย

– ถ ้าเรา Normalized แต่ละ Eigenvector ที่

คํานวณได ้ให ้มี Magnitude เท่ากับหนึ่ง เราจะ

ได ้ Orthogonal Matrix ที่สามารถ Diagonalize

Matrix ที่ต ้องการ

Break

Note: เราเลือก 2 vector ที่

Independent ก ัน โดยการแทน

ค่า β และ γ อะไรก็ได้

End of Chapter 3

• HW 3 Due Next Week

• Chapter 4 Probability Review

Chapter 4 Probability

• Concept and Definition

– Experiments, Event, Outcomes, Sample Space

– Laplace Definition, Definition from Set

• Independent and Mutual Exclusive

– Axioms of Probability, Vein Diagram, Independent

concept, ME Concept

• Conditional Probability and Bayes

• PDF and CDF Concept and Properties

– Continuous and Discrete

Definition

• Outcome/Sample Point

– ผลลัพท์ที่ได ้จากการทดลอง หรือสุ่มตัวอย่าง

• Sample Space

– Set ของผลลัพธ์ทั้งหมด

• Event

– เงื่อนไขของการทดลอง

• กําหนด Event A, ถ ้าทดลอง N ครั้ง และได ้

ผลลัพธ์เป็นไปตามเงื่อนไข A = NA ครั้ง

N

P[ A

A

=

]

P

robabilit y of

E

vent

O

A ccur = lim

N →∞ N

Example: Dice Roll

• Sample Space = {1,2,3,4,5,6}

• P(1)=P(2)= … =P(6)=1/6

– Laplace Definition of Probability

– ถ ้าแต่ละ Member ใน Sample Space มีโอกาส

เกิดเท่าๆกัน

• A = Even

• A = {2,4,6}

• P(A) = |{2,4,6}|/|S| = 3/6 = ½

Mutually Exclusive

• A ME B

• เมื่อเกิด Event A จะเกิด Event B ไม่ได ้

• เมื่อทอยลูเต๋าได ้เลขคู่ จะไม่ใช่เลขคี่

– ดังนั้นถ ้า A = Even, B = Odd

– A ME B

• P(A+B) or P(A union B) = P(A) + P(B)

• เห็นได ้ชัดโดยแสดงด ้วย Vein Diagram

A

B

A

B

Mutually Exclusive

• ME

A

B

• A∩B =∅

• A∪B=A+B

• P(A∪B)=P(A)+P(B)

• Non ME

A

B

• A∩B ≠∅

• A∪B=A+B- A∩B

– Inclusion-Exclusion

Principle

3 AXIOMS OF Probability

• 1. P(S) = 1

– ทุกๆการทดลอง ผลลัพธ์ต ้องอยู่ใน Sample

Space

• 2. 0 =< P(A) =< 1

– ค่าของ Probability ต ้องอยู่ระหว่าง 0 และ 1

• 3. ME: P(A+B)=P(A)+P(B)

– Mutual y Exclusive; Probability ของ Union

ของ Event เท่ากับผลบวกของ Probability

ของแต่ละ Event

Conditional Probability

• Probability ของ Event หนึ่ง เมื่อกําหนดให ้

อีก Event หนึ่งได ้เกิดขึ้น

– Probability จะเพิ่มถ ้าสอง Event เกี่ยวข ้องกัน

– Probability จะไม่เปลี่ยนถ ้าสอง Event ไม่เกี่ยว

กัน

• เราเรียกว่าเป็น Statistical Independent

• นอกจากนี้

[

P A \ B] = [

P AB]/ [

P B],

[

P B \ ]

A = [

P

}

AB / [

P

]

A

[

P A \ B] [

P B] = [

P B \ ]

A

[

P

]

A

[

P A \ B] = [

P B \ ]

A

[

P

]

A / [

P B]

Bayes Rule

E1

E3

E6

A

E

E

9

2

E

4

E7

E

E

8

5

Properties

Example 1

• ในการส่งข ้อมูลแบบ Digital เป็น Frame ขนาด 50 บิต การส่งจะ

สมบูรณ์ได ้ก็ต่อเมื่อทั้ง Frame ไปถึงอย่างถูกต ้อง ถ ้าการเกิด Error ใน

แต่ละบิตของการส่ง(BER = Bit Error Rate) มีโอกาสจะผิดพลาดได ้

เท่ากับ 1/1000 จงหาว่า Frame ที่ส่งจะมี Error เฉลี่ยแล ้วกี่เปอร์เซ็นต์

(FER = Frame Error Rate) ถ ้าการเกิด Error แต่ละ Bit เป็น

Independent

[

P Bit Error] = 10-3

[

P Bit No E

t rror] = 1 -10-3 =

999

.

0

[

P Frame No E

t rror] = P

B

50

[

it No E

t rror] = 0.99950

[

P Frame E

rror] = 1 - 0.99950 = .

0 04879 =

%

879

.

4

Example 2

• จากสถิติ พบว่าโอกาสที่นักศึกษาชายจะสอบผ่านวิชา CPE332 มีค่า

เท่ากับ 0.5 และโอกาสที่นักศึกษาหญิงจะสอบผ่านมีค่าเท่ากับ 0.4 ถ ้า

วิชา CPE332 เทอมนี้มีนักเรียนชาย 40 คนและนักเรียนหญิง 25 คน จง

หาว่าเฉลี่ยแล ้วจะมีนักเรียนตกกี่คน

• สามารถคํานวนโดยใช ้ค่าถ่วงนํ้าหนัก

• ให ้ Sample Space ประกอบด ้วยนักเรียนชาย(M) และนักเรียนหญิง(F)

ดังนั้น S = {M,F}

– สังเกตุว่า Set ทั้งสองเป็น ME คือ Sample Space ถูก Partition

เป็นสอง Partition

– P[M]=40/65 และ P[F]=25/65

• ให ้ Event A เป็นเหตุการณ์ที่นักศึกษาจะสอบผ่าน เราได ้

• P[A\M]=0.5 และ P[A\F]=0.4

• จากสมการการ Partition

n

[

P

]

A = ∑ [

P E ] [

P A \ E ] = [

P M ] [

P A \ M ] + [

P F ] [

P A \ F ]

i

i

i 1

=

40

25

=

⋅ .

0 5 +

⋅ 4

.

0

= 4615

.

0

65

65

[

P Fail C

PE332] = P[AC ] = 1 − [

P

]

A = 1 −

4615

.

0

= 5385

.

0

=

%

85

.

53

Example3

• ต่อจากตัวอย่างที่ 2: ถ ้าเราสุ่มตัวอย่างนักศึกษามาหนึ่งคนที่ลงวิชา CPE332 เมื่อเทอมที่

แล ้ว จากนั้นถามว่านักศึกษาผ่านวิชานี้หรือไม่ นักศึกษาผู ้นั้นตอบว่าสอบผ่านแล ้ว จง

คํานวณว่า

– 1.โอกาสที่นักศึกษาผู ้นั้นจะเป็นผู ้หญิงเท่ากับเท่าไร

– 2. Probability ที่นักศึกษาสอบผ่านและเป็นผู ้หญิงมีเท่าไร

• 1. ต ้องการหา P[F\A]

[

P F \ ]

A = [

P

]

FA / [

P

]

A = [

P A \ F ] [

P F ]/ [

P

]

A

25

1

= 4

.

0 ⋅

= 333

.

0

4 =

%

34

.

33

65 0 461

.

5

• 2. ต ้องการหา P[FA]=P[F∩A]

[

P F \ ]

A = [

P

]

FA / [

P

]

A

[

P

]

FA = [

P F \ ]

A

[

P

]

A = [

P A \ F ] [

P F ]

= .

0 3334 ⋅ 4615

.

0

= .

0 1539 =

%

39

.

15

Random Variables

• เมื่อเรากําหนดค่าเป็นตัวเลขของทุกๆ Sample

Point ใน Sample Space และถ ้าให ้ Variable

แทนผลลัพธ์ที่ได ้จากการทดลอง ดังนั้นผลของการ

ทดลองจะมีค่าเป็นตัวเลข และเราเรียก Variable

นั้นว่าเป็น Random Variable ซึ่งปกติแล ้วเรา

มักจะให ้ Random Variable แทนด ้วยตัว Capital

• ที่สําคัญคือ การกําหนด RV ซึ่งเป็นตัวเลขทําให ้เรา

สามารถนําไปคํานวณต่อทางคณิตศาสตร์ได ้

– Mean

– Variance

– Etc.

Random Variables

• เมื่อผลลัพธ์ของการทดลอง(Sample Space)ได ้ Infinite

Set การกําหนดตัวเลขมักจะเป็นตัวเลขที่ต่อเนื่อง

– เราได ้ Continuous Random Variable

• เมื่อผลลัพธ์เป็น Finite Set การกําหนดจะใช ้ Set ของ

ตัวเลข มักจะเป็น Integer

– เราได ้ Discrete Random Variable

• Probability ที่จะได ้ผลลัพธ์การทดลองหนึ่งๆ คือ

Probability ที่ Random Variable จะมีค่าตามที่เรากําหนด

• เราสามารถแสดงคุณสมบัติของ RV จากการ Plot ค่า

Probability(y-axis) และค่าของ Random Variable(x-

axis)

– Cumulative Distribution Function (CDF)

– Probability Density Function (PDF)

CDF:Cumulative Distribution

Function of RV X

FX(x)

1.0

F (10) = P[ X ≤ 10]

X

x

x = 10

CDF of Normal (Gaussian) Distribution: Continuous

CDF Properties

CDF การทอยเหรียญ: Discrete RV

• ให ้ “หัว”= 0 และ “ก ้อย” = 1

F (x) = P(X ≤ x)

X

1.0

0.5

0

x

0

1

CDF การทอยลูกเต๋า: Discrete RV

• ให ้ ผลเป็นตัวเลขตามหน ้าลูกเต๋า

F (x) = P(X ≤ x)

X

1.0

0.5

0

x

0

1

2

3

4

5

6

7

PDF: Probability Density

Function

f(x)

Area = f(x)dx = 1

x

PDF of Normal(Gaussian) Distribution : Continuous

Properties of PDF

Discrete Version

– ค่าของ Variable ไม่ต่อเนื่อง

• RV X มีค่าเฉพาะที่ X=xi

– F(x) = P(X≤x)

• Function นี้มีความต่อเนื่องด ้านขวามือ

• นิยามสําหรับทุกจุดใน Domain ของ x

• Function เป็นลักษณะขั้นบรรได

• Monotonic Increasing Function จาก 0 ถึง 1

– f(xi) = P(X = xi)

• นิยามเฉพาะจุด ไม่ต่อเนื่อง ค่าเป็นศูนย์ระหว่างนั้น

• บางทีเรียก Probability Mass Function

• ∑f(xi)=1 เสมอ

f(x)

PMF

x

F(x)

CDF

x

Chapter 4 Cont. Next week

– Statistical Average

– Importance PDF

– Joint PDF

– Correlation/Covariance

• Chapter 5: Introduction to Random

Process

– Concept/Definition

– Stationary Concept

– Ergodic Concept

– Autocorrelation and Cross Correlation

CPE 332

Computer Engineering

Mathematics II

Part II, Chapter 4 Cont.: Statistical Average

Chapter 5: Random Process

Today Topics

• HW 3 Due/ HW4 Due Next Week

• PDF

• Expectation

• Joint Density

• Correlation/Covariance

• Random Process

• Stationary

• Ergodic

Random Variables

• Mapping ผลลัพธ์เป็นตัวเลข (One-to-

One) ดังนั้นผลการทดลองสามารถไป

คํานวณต่อได ้

– ที่สําคัญคือค่าเฉลี่ยทางสถิติ

• Mean

• Variance

• Etc.

– และกราฟแสดงคุณสมบัติ Probability

• Cumulative Distribution Function (CDF)

• Probability Density Function (PDF)

– เมื่อผลลัพธ์ของการทดลอง(Sample Space)ได ้ Infinite Set เรา

ได ้ Continuous Random Variable

– เมื่อผลลัพธ์เป็น Finite Set การกําหนดจะใช ้ Set ของตัวเลข

มักจะเป็น Integer เราได้ Discrete Random Variable

CDF:Cumulative Distribution

Function of RV X

FX(x)

1.0

F (10) = P[ X ≤ 10]

X

x

x = 10

CDF of Normal (Gaussian) Distribution: Continuous

CDF Properties

CDF การทอยเหรียญ: Discrete RV

• ให ้ “หัว”= 0 และ “ก ้อย” = 1

F (x) = P(X ≤ x)

X

1.0

0.5

0

x

0

1

PDF: Probability Density

Function

f(x)

Area = f(x)dx = 1

x

PDF of Normal(Gaussian) Distribution : Continuous

Properties of PDF

Discrete Version

– ค่าของ Variable ไม่ต่อเนื่อง

• RV X มีค่าเฉพาะที่ X=xi

– F(x) = P(X≤x)

• Function นี้มีความต่อเนื่องด ้านขวามือ

• นิยามสําหรับทุกจุดใน Domain ของ x

• Function เป็นลักษณะขั้นบรรได

• Monotonic Increasing Function จาก 0 ถึง 1

– f(xi) = P(X = xi)

• นิยามเฉพาะจุด ไม่ต่อเนื่อง ค่าเป็นศูนย์ระหว่างนั้น

• บางทีเรียก Probability Mass Function

• ∑f(xi)=1 เสมอ

f(x)

PMF

x

F(x)

CDF

x

Statistical Average

Discrete V

ersion :

E[ g( X )] = ∑ g( x ) p( x ) i

i

x X

i

Importance Expectation

Discrete V

ersion :

N

1 N

[

E X ] = ∑ x p( x ) = ∑ x p( x ) = ∑ x | p( x ) 1

=

i

i

i

i

x

∀ ∈ X

i=

N

i

i

N

1

i 1

=

i

N

2

2

1

[

E X ] = ∑ x p( x )

2

=

x | p( x ) 1

=

i

i

x

∀ ∈ X

N

i

i

N

i 1

=

i

2

σ = ∑( x − µ )2 p( x ) = [ 2

E X ] − [

E X ]

X

i

X

i

x

∀ ∈ X

i

• Note:

– ในกรณีของ Data ที่ได ้จากการสุ่มตัวอย่าง เราได ้เฉพาะค่า

Estimate ของ Mean และ Mean Square เท่านั้น และ

From n

S

amples : p( x ) = 1 , Each sample h

as the s

ame p

robabilit y

i

N

N

ˆ

µ = X = 1

x

X

N i

i =1

N

Estimate E[ X 2 ] = X 2 = 1

x 2

N i

i =1

2

σ = E[ X 2] − E[ X fr

] om e

stimation i

s bi

ased.

X

We u

se the l

east biased f

ormular :

1

2

2

N

2 

SampleVari ance(Estim ated)

V

:

= s

=

( x

X )

X

N

1

∑ −

N − 1

i

i =1

Multivariate: RV มากกว่า 1 ตัว

• จุดประสงค์ในการศึกษา RV มากกว่าหนึ่งตัว

พร ้อมๆกัน เพื่อจะหาความสัมพันธ์ระหว่าง RV

– แสดงโดยกราฟ Joint PDF

– หรือค่าเฉลี่ยทางสถิติ Correlation และ Covariance

กรณีของ Bivariate

– เป็นการศึกษาความสัมพันธ์ระหว่าง RV สองตัว คือ X

และ Y แต่ละตัวอาจเป็น Discrete หรือ continuous

• เช่นความสัมพันธ์ระหว่างส่วนสูงและอายุ

• หรือความสัมพันธ์ระหว่างคะแนน Midterm และเกรดปลาย

เทอม

• หรือความสัมพันธ์ระหว่าง GPA กับเงินเดือนที่ได ้เมื่อจบ

การศึกษา

Multivariate: RV มากกว่า 1 ตัว

Multivariate: RV มากกว่า 1 ตัว

Joint CDF/Marginal

Density

Joint Moment

Correlation

• G(x,y) = XY

N

M

N

M

Discrete Version : R

= ∑∑ x y p( x , y ) = ∑∑ x y p( x ) p( y ) XY

i

j

i

j

i

j

i

j

i 1

= j 1

=

i 1

= j 1

=

N

M

N

M

1

1

= ∑ x p( x )∑ y p( y ) | independent = ∑ x

y | p( x 1

) =

, p( y

1

) =

i

i

i

i

i

i

i

N

i

M

1

=

1

=

N

1

=

M

i

i

i

i 1

=

Covarience

Correlation and Covarience

การประมาณค่า Rxy และ Cxy จาก

Samples

• เราเก็บตัวอย่างเป็นคู่ (xi,yi) จํานวน N คู่

• Probability ของการได ้แต่ละตัวอย่างเท่ากัน คือ

1/N

• ดังนั้น

N

N

1

R

= E[ XY] =

x y P( x , y )

x y 1

( )

x y

XY

∑∑

=

i

j

i

j

=

i

i

N

i

i

1

N

i

j

i =

i =1

• และ

N

C

= E[( X m ) Y

( − m )] = ∑( x m )( y m

1

)( )

XY

X

Y

i

X

i

Y

N

i=1

N

= 1 ∑( x m )( y m ) = R m m i

X

i

Y

XY

X

Y

N i=1

EX.จงหา Rxy และ Cxy จาก

ตัวอย่างในตารางข ้างล่าง

• ข ้อมูลจากชายไทย อายุระหว่าง 12 – 30 ปี

• X=อายุ และ Y=ส่วนสูง

No.(i) xi

yi No.(i) xi

yi

1

20

178

7

12

168

2

25

176

8

18

156

3

19

163

9

26

174

4

21

184

10

20

171

5

30

180

11

24

182

6

16

165

12

17

179

Scatter Diagram

Calculation Table

PDF ที่สําคัญ

P[X=x]

P[X<=x]

1.0

q

p

x

0

1

Binomial Distribution

n

n

n

!

n

C( ,

n k) = C = C =

=

k

k

 

k k (

! n k)!

b(k;10,0.5)

b(k;10,0.2)

Geometric Distribution

P=0.5

P=0.2

Uniform Distribution

f(x)

1/(b-a)

a

b

x

Gaussian Distribution

Gaussian

Area = 1

Jointly Gaussian: X, Y

P = 0

Volumn = 1

P = 0.5

P = 0.9

P = 0.95

P = -0.99

Exponential Distribution

Exponential Distribution

A = 1

Area = 1

Poisson Distribution

Lambda = 1

Lambda = 3

Lambda = 5

Lambda = 8

Lambda = 12

Lambda = 18

Break

• After Break Chapter 5 Random

Process

Random Process

• เมื่อ Random Variable เป็น Function กับเวลา

– สัญญาณที่มีลักษณะ Random เช่น Noise

– การส่ง Packet ใน Network

– จํานวนรถที่วิ่งบนถนน

– จํานวนลูกค ้าที่เข ้ามาใช ้บริการ

• เหล่านี้ ในแต่ละเวลาหนึ่งๆ ค่าของมันจะมีลักษณะเป็น

Random คือมี PDF ตามแต่ชนิดของมัน

– ที่เวลาต่างกัน PDF อาจจะเปลี่ยน และค่าทางสถิติจะเปลี่ยนตาม

• Random Variable ที่เป็น Function (Random Function

หรือ Analytic Function) กับเวลา เราเรียก Random

Process หรือ Stochastic Process

– จาก RV X จะเป็น X(t)

– Mean และ Variance จะเป็น Function กับเวลาด ้วย

• การวิเคราะห์จะซับซ ้อนกว่า RV

Random Process

(Stochastic Process)

Statistical Average

2

2

2

σ → σ

( t)

X

X ( t )

Sample Function: Ensemble

Average and Time Average

Sample Function and

Ensemble

X ( t)

X ( t )

X ( t )

2

1

ω ( t)

n

Time

ω ( t)

j

Average

1 ∫ ω t() dt

i

T

T

ω ( t)

i

ω ( t)

time

2

ω ( t)

1

t

Ensemble

Ensemble

Average

Average

X ( t),

2

σ ( t)

X

2

X ( t ,

) σ

2

X ( t ,

) σ

1

X ( t )

1

2

X ( t )

2

Statistical Average:

Correlation/Covariance

Stationary RP

2

2

2

σ

( t) → σ

= σ

X ( t )

X ( t )

X

Stationary RP

• Stationary RP

– Strict Sense Stationary (SSS)

• PDF ไม่เปลี่ยน ดังนั้นทุกๆ Moment

(Expectation) ไม่เปลี่ยน

• ปกติ SSS จะหาได ้ยากในทางปฎิบัติ บ่อยครั้งที่เรา

สนใจแค่สอง Moment แรก

– Wide Sense Stationary (WSS)

• เฉพาะ 2 Moment แรกคงที่ แต่ PDF อาจแปรผัน

– Mean

– Variance (Mean Square)

WIDE SENSE

STATIONARY

0

t2-t1 t1

t2

Autocorrelation(WSS)

R

( t , t ) = R ( t t ) = R (τ ) XX

1

2

X

2

1

X

Peak=Power

Even

τ = t t

2

1

Ergodic Random Process

• ถ ้า RP เป็น Ergodic

– Time Average = Ensemble Average

– Mean Ergodic

• หมายถึงค่าเฉลี่ยในทางเวลา 1 ∫ ω t() dt

i เท่ากับ

ค่าเฉลี่ย(mean = E(X(t

T

T

i))) ของ Ensemble

– Correlation Ergodic

• หมายถึงค่า Variance ในทางเวลา 1 ∫ ( ω t() ω t 2

( )) dt

i

i

T

T

เท่ากับค่า Variance ที่ได ้จาก Ensemble 2

σ X( t ) i

– เราใช ้ควบคู่กับ Concept ของ Stationary Random

Process

Discrete-Time Random

Process

Discrete-Time Random

Process (Limited Samples)

Ergodic Discrete RP

• ถ ้า RP เป็น Ergodic ด ้วยค่า Mean (Average

จาก Ensemble = Average จาก Time)

Ergodic Discrete RP

• ถ ้า RP เป็น Ergodic สําหรับ Coorelation

กรณีที่จําก ัด Sequence ความยาว N

จําก ัด Sequence ความยาว N

• ค่า Autocorrelation สําหรับ N Samples

– สมการจะลดรูป เหลือแค่ Sum และเฉลี่ย N

Point แต่จะเกิดการ Biased

N

1

R

( m) = lim

x( n) x( n

m)

XX

+

N →∞ 2 N + 1 n=− N

ความหมายของ x(n)x(n+m)

N

1

R

( m) = lim

x( n) x( n

m)

XX

+

N →∞ 2 N + 1 n=− N

X(n)

X(n+2)

X(n-3)

ความหมายของ x(n)x(n+m), N Points

N

N m 1

R

( m = 1

)

x n x n

m

x n x n

m

XX

∑ ( ) ( + = 1

)

− +

( ) ( +

)

N

N

n=0

n=0

X(n): n=0,1,…,9

X(n+3)

เราไม่ควรเฉลี่ยทั้ง N ตัว (Biased) แต่ควรเฉลี่ย N-m ตัว (non-Biased)

จําก ัด Sequence ความยาว N

• ค่า Autocorrelation สําหรับ N Samples

– สมการจะลดรูป เหลือแค่ Sum และเฉลี่ย N

Point แต่จะเกิดการ Biased เพราะเราเฉลี่ย

น ้อยกว่านั้น

Sequence ความยาว N / RAW

Sequence ท่ัวไป

Example: Rxx

• x=[1 2 3]

Example: Rxx

• x=[1 2 3]

Example: Rxx

• x=[1 2 3]

Example: Rxx

• x=[1 2 3]

Example: Rxx

• x=[1 2 3]

Example: Rxx

• x=[1 2 3]

Example: Rxx

• x=[1 2 3]

Example: Rxx

• x=[1 2 3]

Rxx(m)

Example: Rxy

• x=[1 2 3], y=[4 5 6 7]

Example: Rxy

• x=[1 2 3], y=[4 5 6 7]

Example: Rxy

• x=[1 2 3], y=[4 5 6 7]

Example: Rxy

• x=[1 2 3], y=[4 5 6 7]

Example: Rxy

• x=[1 2 3], y=[4 5 6 7]

Example: Rxy

• x=[1 2 3], y=[4 5 6 7]

Example: Rxy

• x=[1 2 3], y=[4 5 6 7]

Rxy(m)

End of Week 5

• Download HW 4: Probability

– ส่งต ้นชั่วโมงสัปดาห์หน้า

CPE 332

Computer Engineering

Mathematics II

Week 6

Part II, Chapter 5 Random Process

MarKov Process

Chapter 6: Introduction to Queuing

Today Topics

• Random Process

– Stationary: สถิติไม่เปลี่ยน

– Ergodic: Ensemble Average=Time

Average

– Autocorrelation

– Cross Correlation

• ต่อ

• Counting Process

• MarKov Process

Discrete-Time Random

Process

• อาจจะได ้จากการสุ่มตัวอย่าง (Sampling) ของ

Continuous RP

– Uniform Sampling ด ้วย Sampling Period 𝑇𝑇 = 1

𝑓𝑓𝑠𝑠

• เราได ้ … , X −T , 𝑋𝑋 0 , 𝑋𝑋 𝑇𝑇 , 𝑋𝑋 2𝑇𝑇 , 𝑋𝑋 3𝑇𝑇 , …

• ปกติจะละ Sampling Period ไว ้ฐานที่เข ้าใจ เราได ้

… , 𝑋𝑋 −1 , 𝑋𝑋 0 , 𝑋𝑋 1 , 𝑋𝑋 2 , 𝑋𝑋 3 , …

• มักจะเขียนในลักษณะ 𝑋𝑋 𝑛𝑛 ; 𝑛𝑛 = 𝑖𝑖𝑛𝑛𝑖𝑖𝑖𝑖𝑖𝑖𝑖𝑖𝑖𝑖

𝑋𝑋(−8) 𝑋𝑋(−7) 𝑋𝑋(−6) 𝑋𝑋(−5) 𝑋𝑋(−4) 𝑋𝑋(−3) 𝑋𝑋(−2) 𝑋𝑋(−1)

𝑋𝑋(1) 𝑋𝑋(2) 𝑋𝑋(3) 𝑋𝑋(4) 𝑋𝑋(5) 𝑋𝑋(6)

𝑋𝑋(7)

𝑋𝑋(0)

จําก ัด Sequence ความยาว N

• ค่า Autocorrelation สําหรับ N Samples

– สมการจะลดรูป เหลือแค่ Sum และเฉลี่ย N

Point แต่จะเกิดการ Biased เพราะเราเฉลี่ย

น ้อยกว่านั้น

Sequence ท่ัวไป

મ઼��લ઼�

��દ્ધઙ્ખ Correlation

• 1. ด ้วยวิธีกราฟ

– X(n) คูณ X(n+m) คือ X(n) ที่เลื่อนไปซ ้าย m

ตําแหน่ง

• X(n-m) จะเลื่อนไปด ้านขวาแทน

– จากนั้นทําการคูณตัวอย่างที่ตําแหน่งเดียวกัน

และจับผลลัพธ์มาบวกกัน(Summation)

– ถ ้าเป็น Autocorrelation ทําแค่ครึ่งเดียว

เพราะ Rxx(m) เป็น Even Function

• Cross Correlation ต ้องทําทั้งสองด ้าน

મ઼��લ઼�

��દ્ધઙ્ખ Correlation

• 2. ด ้วยการแตก Summation

– 𝑅𝑅

𝑋𝑋𝑋𝑋 𝑚𝑚 = ∑𝑛𝑛=−∞ 𝑥𝑥 𝑛𝑛 𝑥𝑥(𝑛𝑛 + 𝑚𝑚) = 𝑅𝑅𝑋𝑋𝑋𝑋 −𝑚𝑚

– 𝑅𝑅

𝑋𝑋𝑌𝑌 ±𝑚𝑚 = ∑𝑛𝑛=−∞ 𝑥𝑥 𝑛𝑛 𝑦𝑦(𝑛𝑛 ± 𝑚𝑚), not Even

– Index ของ Summation จะสิ้นสุดแค่ Index ของ

x(n) ก็พอ เพราะที่เหลือจะเป็น ศูนย์

– เช่น x(n) = {2, 1, -1 ,3}, y(n) = {1,2,1,-1,-3,4,5}

-1 0 1 2

-3 -2 -1 0 1 2 3

– 𝑅𝑅

2

𝑋𝑋𝑋𝑋 2 = ∑𝑛𝑛=−1 𝑥𝑥 𝑛𝑛 𝑥𝑥 𝑛𝑛 + 2 = 𝑅𝑅𝑋𝑋𝑋𝑋 −2

• = 𝑥𝑥 −1 𝑥𝑥 1 + 𝑥𝑥 0 𝑥𝑥 2 + 𝑥𝑥 1 𝑥𝑥 3 + 𝑥𝑥 2 𝑥𝑥 4

• = 2 −1 + 1 3 + −1 ∙ 0 + 3 ∙ 0 = 1 = 𝑅𝑅𝑋𝑋𝑋𝑋(−2)

– 𝑅𝑅

2

𝑋𝑋𝑌𝑌(−1) = ∑𝑛𝑛=−1 𝑥𝑥 𝑛𝑛 𝑦𝑦 𝑛𝑛 − 1 ≠ 𝑅𝑅𝑋𝑋𝑌𝑌(1)

• = 𝑥𝑥 −1 𝑦𝑦 −2 + 𝑥𝑥 0 𝑦𝑦 −1 + 𝑥𝑥 1 𝑦𝑦 0 + 𝑥𝑥 2 𝑦𝑦 1

• = 2 ∙ 2 + 1 ∙ 1 + −1 ∙ −1 + 3 ∙ −3 = −3

Random Process મ઼�ઙ્ખ �

ઠ્ઠ �

હ્યદ્નદ્બ

• Counting Process

• Birth and Death Process

• Poisson Process

• MarKov Process

• MarKov Chain

Counting Process

N(t)

t

Poisson Process

Poisson Process

• ถ ้าแต่ละเหตุการณ์ที่เกิดขึ้นเป็น Random

และไม่ขึ้นต่อกัน มันจะเป็น Poisson

– Probability ที่จะมี k เหตุการณ์เกิดในช่วงเวลา

t สามารถคํานวณได ้จากสูตร (ดูหน้าถัดไป)

– ระยะเวลาระหว่างสองเหตุการณ์ที่เกิดขึ้น เรียก

Inter-arrival time, 𝜏𝜏, จะมีการกระจายแบบ

Exponential ด ้วยค่าเฉลี่ย 1/λ,

• 𝐹𝐹𝑋𝑋 𝜏𝜏 = 𝑃𝑃 𝑋𝑋 ≤ 𝜏𝜏 = 1−𝜆𝜆𝑒𝑒−𝜆𝜆𝜏𝜏, 𝑓𝑓𝑋𝑋 𝜏𝜏 = 𝜆𝜆𝑖𝑖−𝜆𝜆𝜏𝜏

Poisson Process

Birth and Death Process

State Diagram

• พิจารณาจากระบบ มีทั้ง Birth ด ้วย Birth Rate

λ( t) และ Death ด ้วย Death Rate µ( t)

– เมื่อเราให ้ระบบทํางาน ในระบบจะไม่มีอะไรอยู่ เรา

เรียกว่าอยู่ที่ State 0

– เมื่อมีหนึ่ง Event เข ้ามา หรือ Birth ระบบจะมี Event

เพิ่มขึ้นและจะไปอยู่ที่ State ที่มากกว่าปัจจุบัน “หนึ่ง”

– เมื่อมีหนึ่ง Event จบลง (Death) ระบบจะลด State ลง

หนึ่ง

– ค่า State ของระบบคือจํานวน Event ที่มีอยู่ในระบบ

– การกระโดดไปยัง State ที่สูงกว่า หรือตํ่ากว่า สามารถ

กําหนดด ้วย Probability และเขียนได ้ในลักษณะของ

State Diagram

• ถ ้าระบบไม่มีการจดจํา เราเรียก Diagram นี้เป็น MarKov Model

• ระบบคือ MarKov Process

• การ Transition จาก State หนึ่ง ไปอีก State หนึ่ง กําหนดได ้

โดยค่า Probability ที่คงที่

MarKov Process and Markov Chain

xn-2

x

n-1

x

n

xn+1

xn+2

Markov Process and Markov Chain

n-2

n-1

n

n+1

n+2

MarKov Chain

• สําหรับ MarKov Chain (Discrete State แต่

Time อาจจะเป็น Continuous หรือ Discrete)

– ในขั้นนี้ เราจะเน้นที่ Discrete Time Markov Chain

– สมมุติแกนเวลา ถูกแบ่งเป็น Time Slot ที่เท่ากัน

และเราอยู่ที่ State i ใน Time Slot ปัจจุบัน ดังนั้น

จะมีเหตุการณ์เกิดได ้ดังนี้

• ระบบอยู่ที่ State เดิม ใน Time Slot หน ้า ด ้วยค่า

Probability 𝑃𝑃𝑖𝑖𝑖𝑖

• ระบบมีการเปลี่ยน State ไปยัง State อื่นๆ ด ้วยค่า

Probability 𝑃𝑃𝑖𝑖𝑖𝑖; 𝑗𝑗 = 0,1,2, … , 𝑁𝑁, 𝑗𝑗 ≠ 𝑖𝑖

Note: Continuous Time MarKov Chain สามารถ Model จาก Discrete Time

MarKov Chain โดยให้ Limit ระยะห่างของ Time Slot เข้าสู่ศูนย์

Discrete Time Markov Chain

ในรูปไม่ได ้แสดง Transition หมดทุกเส ้น

𝑃𝑃0,4

𝑃𝑃0,3

𝑃𝑃0,2

𝑃𝑃0,0

𝑃𝑃1,1

𝑃𝑃2,2

𝑃𝑃3,3

𝑃𝑃4,4

𝑃𝑃

𝑃𝑃

𝑃𝑃

𝑃𝑃

0

0,1

3,4

1

1,2

2

2,3

3

4

𝑃𝑃1,0

𝑃𝑃2,1

𝑃𝑃3,2

𝑃𝑃4,3

𝑃𝑃2,0

𝑃𝑃4,0

𝑃𝑃3,0

Discrete Time Markov Chain

ถ้าระบบมีได้ถึง State N, (State 0,1,2,…,N),

Transition Matrix จะมีขนาด (N+1)คูณ(N+1)

Discrete Time Markov Chain

Notations: દ્ભ�દ્મ�� Discrete Time MarKov

Chain

• State Probability

– คือ Probability ที่จะพบว่าระบบอยู่ที่ State ใดๆ

𝑝𝑝𝑖𝑖 = 𝑃𝑃𝑖𝑖𝑃𝑃𝑃𝑃𝑃𝑃𝑃𝑃𝑖𝑖𝑃𝑃𝑖𝑖𝑖𝑖𝑦𝑦 𝑆𝑆𝑖𝑖𝑃𝑃𝑖𝑖𝑖𝑖 = 𝑖𝑖, 𝑖𝑖 = 0,1, … , 𝑁𝑁

– ∑𝑁𝑁𝑥𝑥=0 𝑝𝑝𝑥𝑥 = 1

• Transition Probability

– คือ Probability ที่ระบบจะกระโดดจาก State

หนึ่งใน Time Slot ปัจจุบัน ไปยังอีก State หนึ่ง

ใน Time Slot หน ้า อย่าลืมว่าระบบไม่มีการจํา

– 𝑃𝑃𝑖𝑖𝑖𝑖 = 𝑃𝑃 𝑋𝑋𝑛𝑛 = 𝑖𝑖, 𝑋𝑋𝑛𝑛+1 = 𝑗𝑗 ; 𝑖𝑖, 𝑗𝑗 = 0,1, … , 𝑁𝑁

– ∑𝑁𝑁𝑖𝑖=0 𝑃𝑃𝑖𝑖𝑖𝑖 = 1

• ระบบอยู่ที่ Equilibrium ค่าเหล่านี้จะไม่เปลี่ยน

Discrete Time Markov Chain

Discrete Time Markov Chain

Discrete Time Markov Chain

Discrete Time Markov Chain

Discrete Time Markov Chain

สรุป MarKov Chain

• สถานะของระบบ ดูได ้จากจํานวณ Event ที่อยู่ในระบบ

เรียก State ของระบบ

– ถ ้า State เป็น Discrete เราได ้ MarKov Chain

• มีค่า Probability สองชุดที่อธิบายการทํางานของระบบ

– State Probability: Probability ที่ระบบจะอยู่ที่ State ใด

State หนึ่ง

• ผลรวมของ State Probability จะต ้องเท่ากับ 1

– Transition Probability: Probability ที่ระบบจะมีการเปลี่ยน

State

• อธิบายจาก Transition Matrix

• ผลรวมของ Transition Probability แต่ละแถว จะต ้องเท่ากับ 1

• ถ ้าระบบอยู่ที่ Equilibrium ค่า Probability ของ State

จะไม่เปลี่ยน และสามารถอธิบายได ้ด ้วย Global Balance

Equation

• MarKov Chain ที่เราสนใจคือ Irreducible และ

Aperiodic

Detailed Balance Equation:

Simple MarKov Chain

• Detailed Balance Equation

ระบบจะอยู่ที่ State เด ิม หร ือกระโดดไปย ัง State ข ้างเค ียงเท่าน ั�น

Markov Chain(Detailed Bal Eq)

• Detailed Balance Equation

• Transition Matrix ของ Simple Markov

Chain จะมีลักษณะเป็น Tridiagonal

P

P

0

00

01

P

P

P

10

11

12

P = 

P

P

P

21

22

23

Pn− ,1 n

 0

P

P

n, n−1

nn

Example

• MarKov Chain แสดงด ้วย Markov Model

ดังรูปข ้างล่าง

– 1. จงหา Transition Matrix

– 2. จาก Transition Matrix จงคํานวณหา State

Probability

0.2

0.3

0.1

0.2

0

1

2

3

0.1

0.1

0.2

0.5

0.8

0.6

0.6

0.1

0.2

Example

– 1. จงหา Transition Matrix

0.2

0.3

0.1

0.2

0

1

2

3

0.1

0.1

0.2

0.5

0.8

0.6

0.6

0.1

0.2

P

P

P

P

00

01

02

03 

P

P

P

P

10

11

12

13 

P =  P P P P

 20

21

22

23 

P

P

P

P

30

31

32

33 

Example

– 1. จงหา Transition Matrix

0.2

0.3

0.1

0.2

0

1

2

3

0.1

0.1

0.2

0.5

0.8

0.6

0.6

0.1

0.2

P

P

P

P

00

01

02

03 

 5

.

0

3

.

0

0

0 2

. 

 

P

P

P

P

10

11

12

13 

P =

=  1

.

0

0 8

.

0 1

.

0 

P

P

P

P

 1

.

0

.

0 1

6

.

0

0 2

. 

 20

21

22

23  

P

P

P

P

30

31

32

33 

 0

2

.

0

2

.

0

0 6

. 

Example

– 2. จงคํานวณหา State Probability

0.2

 5

.

0

0 3

.

0

.

0 2

0.3

0.1

0.2

0

1

2

3

P = 0 1

.

8

.

0

0 1

.

0 

0.1

0.1

0.2

0 1

.

1

.

0

0 6

.

.

0 2

0.5

0.8

0.6

0.6

0.1

0.2

 0

2

.

0

2

.

0

.

0 6

Example

– 2. จงคํานวณหา State Probability

0.2

 5

.

0

0 3

.

0

.

0 2

0.3

0.1

0.2

0

1

2

3

P = 0 1

.

8

.

0

0 1

.

0 

0.1

0.1

0.2

0 1

.

1

.

0

0 6

.

.

0 2

0.5

0.8

0.6

0.6

0.1

0.2

 0

2

.

0

2

.

0

.

0 6

p + p + p + p = 1 0

.

0

1

2

3

Global B

alanced Eq

uation :

p

P

p P

j

∑∞ =

ji

∑∞ i ij

i=0, ij

i=0, ij

Example

– 2. จงคํานวณหา State Probability

0.2

 5

.

0

0 3

.

0

.

0 2

0.3

0.1

0.2

0

1

2

3

P = 0 1

.

8

.

0

0 1

.

0 

0.1

0.1

0.2

0 1

.

1

.

0

0 6

.

.

0 2

0.5

0.8

0.6

0.6

0.1

0.2

 0

2

.

0

2

.

0

.

0 6

p + p + p + p =

0

.

1

0

1

2

3

Global B

alanced Eq

uation :

p

P

p P

j

∑∞ =

ji

∑∞ i ij

i=0, ij

i=0, ij

p [ P + P + P ] = p P + p P + p P

0

01

02

03

1 10

2

20

3 30

Example

– 2. จงคํานวณหา State Probability

0.2

 5

.

0

0 3

.

0

.

0 2

0.3

0.1

0.2

0

1

2

3

P = 0 1

.

8

.

0

0 1

.

0 

0.1

0.1

0.2

0 1

.

1

.

0

0 6

.

.

0 2

0.5

0.8

0.6

0.6

0.1

0.2

 0

2

.

0

2

.

0

.

0 6

p + p + p + p =

0

.

1

0

1

2

3

Global B

alanced Eq

uation : p [ P + P + P ] = p P + p P + p P

0

01

02

03

1 10

2

20

3

30

p

P

p P

p [ P + P + P ] = p P + p P + p P

j

∑∞ =

ji

∑∞ i ij

1

10

12

13

0

01

2

21

3

31

i=0, ij

i=0, ij

p [ P + P + P ] = p P + p P + p P

2

20

21

23

0

02

1 12

3

32

p [ P + P + P ] = p P + p P + p P

3

30

31

32

0

03

1 13

2

23

Example

– 2. จงคํานวณหา State Probability

0.2

 5

.

0

0 3

.

0

.

0 2

0.3

0.1

0.2

0

1

2

3

P = 0 1

.

8

.

0

0 1

.

0 

0.1

0.1

0.2

0 1

.

1

.

0

0 6

.

.

0 2

0.5

0.8

0.6

0.6

0.1

0.2

 0

2

.

0

2

.

0

.

0 6

p + p + p + p =

0

.

1

0

1

2

3

Global B

alanced Eq

uation :

p [ 3

.

0

+ 0 +

]

2

.

0

= 1

.

0 p +

1

.

0 p + .

0 p

0

1

2

3

p

P

p P

j

∑∞ =

ji

∑∞ i ij

p [ 1

.

0 +

1

.

0 + ]

0 =

3

.

0 p +

1

.

0 p +

2

.

0 p

i=0, ij

i=0, ij

1

0

2

3

p [ 1

.

0 +

1

.

0 +

]

2

.

0

= .

0 p + .

0 1 p + .

0 2 p

2

0

1

3

p [0 +

2

.

0

+

]

2

.

0

= 2

.

0 p + .

0 p +

2

.

0 p

3

0

1

2

Example

– 2. จงคํานวณหา State Probability

0.2

 5

.

0

0 3

.

0

.

0 2

0.3

0.1

0.2

0

1

2

3

P = 0 1

.

8

.

0

0 1

.

0 

0.1

0.1

0.2

0 1

.

1

.

0

0 6

.

.

0 2

0.5

0.8

0.6

0.6

0.1

0.2

 0

2

.

0

2

.

0

.

0 6

p + p + p + p =

0

.

1

0

1

2

3

Global B

alanced Eq

uation :

− 5

.

0 p + .

0 1 p + 0.1 p

= 0

0

1

2

p

P

p P

j

∑∞ =

ji

∑∞ i ij

.

0 3 p

2

.

0 p +

1

.

0 p +

2

.

0 p = 0

i=0, ij

i=0, ij

0

1

2

3

1

.

0 p

4

.

0 p +

2

.

0 p = 0

1

2

3

2

.

0 p

+ 0 2

. p

4

.

0 p = 0

0

2

3

Homogeneous Linear Eq. ยังแก้สมการไม่ได้

Example

– 2. จงคํานวณหา State Probability

0.2

 5

.

0

0 3

.

0

.

0 2

0.3

0.1

0.2

0

1

2

3

P = 0 1

.

8

.

0

0 1

.

0 

0.1

0.1

0.2

0 1

.

1

.

0

0 6

.

.

0 2

0.5

0.8

0.6

0.6

0.1

0.2

 0

2

.

0

2

.

0

.

0 6

− 5

.

0 p + .

0 1 p +

1

.

0 p

= 0

0

1

2

3

.

0 p − 0.2 p + .

0 1 p + .

0 2 p = 0

0

1

2

3

0 1

. p

4

.

0 p + 0.2 p = 0

1

2

3

p

+ p

+ p

+ p = .

1 0

0

1

2

3

ต ้องนําสมการที่ 4 มาใส่เสมอ + อีก 3 สมการอะไรก็ได ้

สมการ Non Homogeneous System of Linear Equation

แก ้ได ้โดยวิธี Elimination จะได ้ Unique Solution

Example

– 2. จงคํานวณหา State Probability

0.2

 5

.

0

0 3

.

0

.

0 2

0.3

0.1

0.2

0

1

2

3

P = 0 1

.

8

.

0

0 1

.

0 

0.1

0.1

0.2

0 1

.

1

.

0

0 6

.

.

0 2

0.5

0.8

0.6

0.6

0.1

0.2

 0

2

.

0

2

.

0

.

0 6

− 0 5

. p +

1

.

0 p + 0 1

. p

= 0

p

0

1

2

− 5

.

0

1

.

0

1

.

0

0  0 

0

   

.

0 3 p − 0 2

. p + 0 1

. p + .

0 2 p = 0

p

0

1

2

3

⇒  3

.

0

− 2

.

0

0 1

.

.

0 2 1 = 0

1

.

0 p − .

0 4 p + 0 2

. p = 0

 0

1

.

0

− 4

.

0

.

0 2 p

0

1

2

3

 2   

p

+ p

+ p

+ p = 1 0

.

p

0

1

2

3

 1

1

1

1  3  1

Non Homogeneous System of Linear Equation

แก ้ได ้โดยวิธี Elimination จะได ้ Unique Solution

Algorithm ในการแก ้สมการ จะกล่าวภายหล ัง

Example

– 2. จงคํานวณหา State Probability

0.2

 5

.

0

0 3

.

0

.

0 2

0.3

0.1

0.2

0

1

2

3

P = 0 1

.

8

.

0

0 1

.

0 

0.1

0.1

0.2

0 1

.

1

.

0

0 6

.

.

0 2

0.5

0.8

0.6

0.6

0.1

0.2

 0

2

.

0

2

.

0

.

0 6

− .

0 5 p +

1

.

0 p +

1

.

0 p

= 0

)

1

(

0

1

2

( )

4 − 5× ( )

2 : − .

0 5 p + 2 p + 0 5

. p = 1

)

5

(

.

0 3 p − 2

.

0 p +

1

.

0 p + .

0 2 p = 0

(2)

0

1

2

0

1

2

3

( )

4 − 5×

)

3

(

:

p + 0.5 p + 3 p = 1 ( )

6

1

.

0 p − 4

.

0 p +

2

.

0 p = 0

)

3

(

0

1

2

1

2

3

p

+ p

+ p

+ p = 0

.

1

( )

4

จากสมการที่ (1),(5),(6) เราได ้

0

1

2

3

− 0.5 p + 0.1 p + 0.1 p = 0

)

1

(

0

1

2

− 0.5 p + 2 p

+ 0.5 p = 1

)

5

(

0

1

2

p + 0.5 p + 3 p = 1 (6)

0

1

2

Example

– 2. จงคํานวณหา State Probability

0.2

 5

.

0

0 3

.

0

.

0 2

0.3

0.1

0.2

0

1

2

3

P = 0 1

.

8

.

0

0 1

.

0 

0.1

0.1

0.2

0 1

.

1

.

0

0 6

.

.

0 2

0.5

0.8

0.6

0.6

0.1

0.2

 0

2

.

0

2

.

0

.

0 6

− 0 5

. p + .

0 1 p +

1

.

0 p = 0

)

1

(

0

1

2

+ ×

p +

p =

(6) 2

)

1

( :

.

0 7

3 2

.

1 (7)

5

.

0 p + 2 p

+ 5

.

0 p = 1

)

5

(

1

2

0

1

2

(6) + 2 ×

)

5

(

:

.

4 5 p +

4 p = 3

)

8

(

p + .

0 5 p + 3 p = 1 ( )

6

1

2

0

1

2

5

.

4

4 5

.

4 5

.

)

8

(

(7) : (4 −

3 )

2

.

p = 3 −

0 7

.

0 7

.

2

0 7

.

− .

3 4285714

p =

= .

0 2069

2

− .

16 57142857

Example

– 2. จงคํานวณหา State Probability

0.2

 5

.

0

0 3

.

0

.

0 2

0.3

0.1

0.2

0

1

2

3

P = 0 1

.

8

.

0

0 1

.

0 

0.1

0.1

0.2

0 1

.

1

.

0

0 6

.

.

0 2

0.5

0.8

0.6

0.6

0.1

0.2

 0

2

.

0

2

.

0

.

0 6

p = 0 2069

.

2

1− .

3 2 p

From (7)

:

2

p =

= 0.4828

1

.

0 7

From )

1

( : p = 0.1379

0

From ( )

4 : p = .

0 1724

3

สมการ System of Linear Equation แก ้ได ้โดยวิธี Elimination

p =

137

.

0

,

9 p = 0 482

.

,

8 p = 0 206

.

,

9 p = 0 172

.

4

0

1

2

3

Example: Simple MarKov

• MarKov Chain แสดงด ้วย Markov Model

ดังรูปข ้างล่าง

– 1. จงหา Transition Matrix

– 2. จาก Transition Matrix จงคํานวณหา

State Probability

0.5

0.1

0.3

0.2

0

1

2

3

4

0.1

0.2

0.2

0.6

0.5

0.8

0.5

0.6

0.4

Example: Simple MarKov

• Transition Matrix

0.5

0.1

0.3

0.2

0

1

2

3

4

0.1

0.2

0.2

0.6

0.5

0.8

0.5

0.6

0.4

P

P

P

P

P

0,0

0 1

,

0,2

0,3

0,4 

0 5

.

.

0 5

0

0

0 

 

P

P

P

P

P

,

1 0

1

,

1

,

1 2

,

1 3

,

1 4 

 .

0 1

.

0 8

1

.

0

0

0 

P ==  P

P

P

P

P  = 

2,0

2 1

,

2,2

2,3

2,4

0

0 2

.

0 5

.

.

0 3

0 

 

P

P

P

P

P

3,0

3 1

,

3,2

3,3

3,4 

 0

0

2

.

0

6

.

0

0 2

. 

P

P

P

P

P

 4,0

4 1

,

4,2

4,3

4,4 

 0

0

0

.

0 6

.

0 4

Example: Simple MarKov

• State Probability

0.5

0.1

0.3

0.2

0

1

2

3

4

0.1

0.2

0.2

0.6

0.5

0.8

0.5

0.6

0.4

P

P

P

P

P

0,0

0 1

,

0,2

0,3

0,4 

 5

.

0

5

.

0

0

0

0 

 

P

P

P

P

P

,

1 0

1

,

1

,

1 2

,

1 3

,

1 4 

 1

.

0

0 8

.

.

0 1

0

0 

P ==  P

P

P

P

P  = 

2,0

2 1

,

2,2

2,3

2,4

0

2

.

0

5

.

0

.

0 3

0 

 

P

P

P

P

P

3,0

3 1

,

3,2

3,3

3,4 

 0

0

.

0 2

6

.

0

.

0 2

P

P

P

P

P

4,0

4 1

,

4,2

4,3

4,4

0

0

0

.

0 6

4

.

0 

 

Detailed B

alance E

quation : p P = p P

i

i, j

j

j , i

p P = p P

.

0 5 p = .

0 1 p o

r p = 5 p

0 0 1

,

1

,

1 0

0

1

1

0

p P = p P

1

.

0 p = .

0 2 p o

r p = 0 5

. p = .

2 5 p

1

,

1 2

2

2 1

,

1

2

2

1

0

p P

= p P

3

.

0 p = .

0 2 p o

r p = .

1 5 p = .

3 75 p

2

2,3

3 3,2

2

3

3

2

0

1

p P

= p P

2

.

0 p = 0 6

. p o

r p =

p = .

1 25 p

3 3,4

4

4,3

3

4

4

3

3

0

Example: Simple MarKov

• State Probability

Detailed B

alance E

quation : p P = p P

i

i, j

j

j , i

p P = p P

.

0 5 p = .

0 1 p o

r p = 5 p

0 0 1

,

1

,

1 0

0

1

1

0

p P = p P

1

.

0 p = .

0 2 p o

r p = 0 5

. p = .

2 5 p

1

,

1 2

2

2 1

,

1

2

2

1

0

p P

= p P

3

.

0 p = .

0 2 p o

r p = .

1 5 p = .

3 75 p

2

2,3

3 3,2

2

3

3

2

0

1

p P

= p P

2

.

0 p = 0 6

. p o

r p =

p = .

1 25 p

3 3,4

4

4,3

3

4

4

3

3

0

From

p + p + p + p + p = 1

0

1

2

3

4

p + 5 p +

5

.

2 p + 3 75

.

p + .

1 25 p = 1

0

0

0

0

0

p 1

[ + 5 + 2 5

. +

75

.

3

+ .

1

]

25 = 1

0

p = 0 0740

.

74

then

0

p = 0 3703

.

,

7 p =

1851

.

0

,

85 p = .

0 2777

,

78 p = .

0 092593

1

2

3

4

Example: Simple MarKov

• State Probability

Detailed B

alance E

quation : p P = p P

i

i, j

j

j , i

p P = p P

.

0 5 p = .

0 1 p o

r p = 5 p

0 0 1

,

1

,

1 0

0

1

1

0

p P = p P

1

.

0 p = .

0 2 p o

r p = 0 5

. p = .

2 5 p

1

,

1 2

2

2 1

,

1

2

2

1

0

p P

= p P

3

.

0 p = .

0 2 p o

r p = .

1 5 p = .

3 75 p

2

2,3

3 3,2

2

3

3

2

0

1

p P

= p P

2

.

0 p = 0 6

. p o

r p =

p = .

1 25 p

3 3,4

4

4,3

3

4

4

3

3

0

From

p + p + p + p + p = 1

ใช้วิธีของ Matrix

0

1

2

3

4

0 5

. p

1

.

0 p = 0

p

0

1

0 5

.

− 1

.

0

0

0

0  0 

0

   

1

.

0 p − 0 2

. p = 0

p

1

2

 0

1

.

0

− 2

.

0

0

0  1 0

0 3

. p − 0 2

. p = 0

⇒  0

0

3

.

0

− .

0 2

0  p  = 0

2

3

 2  

.

0 2 p

6

.

0 p = 0

p

3

4

 0

0

0

2

.

0

− 6

.

0  3 0

p + p + p + p + p = 1

 1

1

1

1

1  p

 

0

1

2

3

4

 4 1

End of Chapter 5

• Chapter 6

– Introduction to Queuing Theory

• Homework Chapter 5 Download

– เน้นที่ Correlation ของ Stationary RP และ

การใช ้ Global/Detailed Balance Equation

ใน MarKov Chain

CPE 332

Computer Engineering

Mathematics II

Week 6

Part II, Chapter 6

Queuing System

Topics

• Birth and Death Process

• Unlimited Server

• N Servers

• Single Server, M/M/1

• Kendal Notation

• Applications

System

• ระบบประกอบด ้วย Input และ Output

• พิจารณาระบบที่มีการให ้บริการ(Service) แก่

ลูกค ้า (Customer)

– ลูกค ้าเข ้ามาในระบบเพื่อขอรับบริการ (Input)

– ระบบมี Resource ที่จํากัดในการให ้บริการ

– ลูกค ้า เมื่อได ้รับบริการแล ้ว ออกไปจากระบบ

(Output)

• ระบบขายของหน้าร ้าน, ระบบหน้าธนาคาร, ระบบ

การจราจร, ระบบ Operating System ใน

คอมพิวเตอร์, สถานีนํ้ามัน/แก๊ส, คิวจ่ายของ/

อาหาร, ระบบสื่อสารข ้อมูล, ระบบโทรศัพท์ และ

อื่นๆอีกมาก

Queuing System

Arrival Rate = λ

Departure Rate = µ

Customer

System

Customer

1. อ ัตร าการ เข ้ามาของลูกค ้า ค ือ Arrival Rate

2. อ ัตร าการ ออกไปของลูกค ้าเม ื่อได ้ร ับบร ิการ เสร็จค ือ Departure Rate

3. State ของร ะบบค ือจํานวนลูกค ้าที่อย ู่ในร ะบบ ที่ร อบร ิการ

หร ือกําล ังถูกบร ิการ

4. ถ ้าระบบไม่ม ีการจดจํา การบร ิการลูกค ้าแต่ละรายเหม ือนก ัน

และไม่ข ึ�นก ับอด ีต เราสามารถใช ้รูปแบบ MarKov Chain อธิบายระบบ

Queuing System

Birth Rate

Death Rate

Arrival Rate = λ

Departure Rate = µ

Customer

System

Customer

5. ถ ้าเราแบ่งช ่วงเวลาการเข ้ามาของลูกค ้าเป็นช ่วง Time Slot

และบ ันทึกเหตุการณ ์ที่เก ิดข ึ�นในแต่ละ Time Slot

a) ม ีลูกค ้าใหม ่เข ้ามา เพ ื่อขอใช ้บร ิการ ( State ของร ะบบจะเพ ิ่ม) หร ือ b) ม ีลูกค ้าที่ได ้ร ับบร ิการแล้วออกจากระบบไป ( State จะลด) หร ือ c) ไม่ม ีลูกค ้าใหม่ และไม่ม ีลูกค ้าที่ให้บร ิการเสร็จ (State คงเด ิม) 6. จาก Model ข ้อ 5 เราจะได ้ Discrete Time Markov Chain

7. ถ ้าช ่วงเวลาของ Time Slot ส ั�นมากจนกร ะท ่ังลูกค ้าที่เข ้ามา หร ือออกไป

ในช ่วงหน ึ่ง Time Slot สามารถม ีได ้แค่คนเด ียว

เราจะสามารถ Model ระบบได ้เป็น Simple MarKov Model

Queuing System

ถ ้า λ < µ ระบบจะสามารถ

เข ้าสู่ Equilibrium ได ้

Birth Rate

; Server Utilization Death Rate

𝜌𝜌 = 𝜆𝜆𝜇𝜇

Arrival Rate = λ

Departure Rate = µ

Customer

System

Customer

Queuing System

Simple MarKov Model

Birth Rate

Death Rate

λ < µ

Arrival Rate = λ

Departure Rate = µ

Customer

System

Customer

Detailed Balance Equation: 𝒑𝒑𝒊𝒊𝑷𝑷𝒊𝒊𝒊𝒊 = 𝒑𝒑𝒊𝒊𝑷𝑷𝒊𝒊𝒊𝒊

สําหร ับสอง State i, j ที่อยู่ต ิดก ันใดๆ

Queuing System

Simple MarKov Model

𝒑𝒑

Birth Rate

𝒊𝒊𝑷𝑷𝒊𝒊𝒊𝒊 = 𝒑𝒑𝒊𝒊𝑷𝑷𝒊𝒊𝒊𝒊

Death Rate

λ < µ

Arrival Rate = λ

Departure Rate = µ

Customer

System

Customer

ในที่นี�เราจะพ ิจารณ ากรณ ี

1.

ที่ลูกค ้าแต่ละรายเข ้ามาในระบบแบบ Random

และเป็น Poisson Process

2.

เวลาในการให้บร ิการของลูกค ้าแต่ละราย เป็น Random

ม ีการ กร ะจายแบบ Exponential

Queuing System

Simple MarKov Model

𝒑𝒑

Birth Rate

𝒊𝒊𝑷𝑷𝒊𝒊𝒊𝒊 = 𝒑𝒑𝒊𝒊𝑷𝑷𝒊𝒊𝒊𝒊

Death Rate

λ < µ

Arrival Rate = λ

Departure Rate = µ

Customer

System

Customer

Arrival: 1. Probability ที่จะม ีลูกค ้า k คนเข ้ามาในช ่วง T ว ินาที

(𝝀𝝀𝝀𝝀)𝒌𝒌𝒆𝒆𝝀𝝀𝝀𝝀

𝑷𝑷 𝒌𝒌 = 𝑷𝑷 𝑿𝑿 = 𝒌𝒌 =

𝒌𝒌!

; 𝒌𝒌 = 𝟎𝟎, 𝟏𝟏, 𝟐𝟐, …

เม ื่อ 𝝀𝝀 เป็นค ่าเฉลี่ยจํานวนลูกค ้าที่เข ้ามา ต ่อว ินาที

2. ค่า Inter-arrival Time, 𝝉𝝉 จะมีการกระจายแบบ

Exponential ด ้วยค ่าเฉล ี่ย 𝟏𝟏 𝝀𝝀

𝑷𝑷 𝝀𝝀 ≤ 𝝉𝝉 = 𝟏𝟏 − 𝝀𝝀𝒆𝒆−𝝀𝝀𝝉𝝉

Queuing System

Simple MarKov Model

𝒑𝒑

Birth Rate

𝒊𝒊𝑷𝑷𝒊𝒊𝒊𝒊 = 𝒑𝒑𝒊𝒊𝑷𝑷𝒊𝒊𝒊𝒊

Death Rate

λ < µ

Arrival Rate = λ

Departure Rate = µ

Customer

System

Customer

Departure: 1. เวลาเฉล ี่ยที่ลูกค ้าใช ้บร ิการ = 𝝀𝝀𝒔𝒔 (Service Time)

จะม ีการกระจายแบบ Exponential

Probability ที่ลูกค ้าจะใช ้เวลาบร ิการ น ้อยกว่าหร ือ

เท่าก ับ 𝒕𝒕 𝑷𝑷 𝝀𝝀 ≤ 𝒕𝒕 = 𝟏𝟏 − 𝟏𝟏 𝒆𝒆−𝒕𝒕 𝝀𝝀�𝒔𝒔

𝝀𝝀𝒔𝒔

2. Departure Rate ค ืออ ัตร าที่ลูกค ้าออกจากร ะบบ

เม ื่อได ้ร ับบร ิการเสร็จ(อ ัตราการให้บร ิการแก่ลูกค ้า)

𝝁𝝁 = 𝟏𝟏

𝝀𝝀𝒔𝒔

Queuing System Case 1:

Unlimited Server; No Queue

Arrival Rate = λ

Departure Rate = µ

Customer

System

Customer

1.

ลูกค ้าที่เข ้ามา เป็น Random ด ้วยอัตราเฉล่ีย 𝜆𝜆 และมีการกระจายแบบ Poisson 2.

สมมุติว่า Customer แต่ละคนที่เข ้ามาได ้รับการ Service จากระบบทันที (ระบบมี

Server จํานวนไม่จํากัด และรับ Customer ได ้ไม่จํากัด)

3.

เวลาที่ใช ้ในการ Service เป็น เป็น Exponential ด ้วยเวลาเฉล่ีย 𝑇𝑇𝑠𝑠

4.

ระบบสามารถรับ Customer ได ้ไม่จํากัด

5.

ลูกค ้าเข ้ามาได ้ทีละคน และออกทีละคน

6.

ระบบนี้เรียก M/M/∞

Queuing System Case 1:

Unlimited Server; No Queue

λ < µ

k

−λ

λ

x

e

ρ ρ

e

[

P k] =

λ

p =

; ρ =

t / Ts

1

[

P T

t] = e

; T =

k!

x

µ

s

µ

!

x

Arrival Rate = λ

Departure Rate = µ

Customer

System

Customer

1.

สมมุติว่า Customer แต่ละคนที่เข ้ามาเป็น Poisson และได ้รับการ Service จากระบบทันที

2.

เวลาที่ใช ้ในการ Service เป็น Random สมมุติว่าเป็น Exponential ด ้วยเวลาเฉลี่ย T

3.

ระบบสามารถรับ Customer ได ้ไม่จํากัด แต่เข ้ามาได ้ทีละคน และออกทีละคน

4.

ระบบนี้เรียก M/M/∞ แสดงได ้ด ้วย Simple Markov Model

5.

เราสามารถพิสูจน์ได ้ว่าค่า State Probability จะมีการกระจายแบบ Poisson

p P = p P

i

ij

j

ji

0

1

2

i

j

Queuing System Case 2: Lost System

Limited Server=N; No Queue

λ < µ

x

ρ / x!

k

−λ

λ

p =

e

x

N

k

ρ

[

P k] =

t / Ts

1

[

P T

t] = e

; T =

k!

s

µ

!

0 k

k =

Arrival Rate = λ

Departure Rate = µ

Customer

System, N Server

Customer

1.

สมมุติว่า Customer แต่ละคนที่เข ้ามาเป็น Poisson และได ้รับการ Service จากระบบทันที

2.

เวลาที่ใช ้ในการ Service เป็น Random สมมุติว่าเป็น Exponential ด ้วยเวลาเฉลี่ย T

3.

ระบบสามารถรับ Customer ได ้ N ถ ้าทุก Server เต็ม จะรับ Customer ใหม่ไม่ได ้ (Lost) 4.

ระบบนี้เรียก M/M/N/N แสดงได ้ด ้วย N-state Simple Markov Model

5.

State Probability จะมีการกระจายแบบ First Erlang (Erlang B) Distribution

p P = p P

i

ij

j

ji

0

1

2

i

j

N

Queuing System Case 3: Delay System

Limited Server=N; With Unlimited Queue

x

ρ

λ < µ

p ; 0 ≤ x N

1

0

N

1



k

!

x

Nρ

ρ 

k

−λ

λ

p =

p

x

x

=

+

N

N

e

; 0

ρ

[

P k] =

N

N!( N − ρ)

k

k =

!

t / T

0

s

1

  p ; x N

[

P T

t] = e

; T =

k!

s

0

µ

 N!  N

Arrival Rate = λ

Departure Rate = µ

Customer

System, N Server

Customer

1.

สมมุติว่า Customer แต่ละคนที่เข ้ามาเป็น Poisson และได ้รับการ Service จากระบบทันที

2.

เวลาที่ใช ้ในการ Service เป็น Random สมมุติว่าเป็น Exponential ด ้วยเวลาเฉลี่ย T

3.

ระบบสามารถรับ Customer ได ้ไม่จํากัด แต่จะ Service ได ้สูงสุด N พร ้อมๆกัน

4.

ถ ้าทุก Server เต็ม Customer ใหม่จะต ้องรอใน Queue ในกรณีนี้จะเกิด Queuing Delay 5.

ระบบนี้เรียก M/M/N หรือ M/M/N/∞ แสดงได ้ด ้วย Simple Markov Model

6.

State Probability จะมีการกระจายแบบ Second Erlang (Erlang C) Distribution

0

1

2

N

N+1

p P = p P

i

ij

j

ji

Queuing System Case 3: Delay System

Server=1; With Unlimited Queue; M/M/1

λ < µ 1

 ρ

−

p =

+1

=1− ρ

0

k

−λ

λ

 − ρ

e

1

(

)

[

P k] =

t / Ts

1

x

[

P T

t] = e

; T =

k!

s

p = ρ 1

( − ρ)

µ

x

Arrival Rate = λ

Departure Rate = µ

Customer

System, 1 Server

Customer

1.

สมมุติว่า Customer แต่ละคนที่เข ้ามาเป็น Poisson และได ้รับการ Service จากระบบทันที

2.

เวลาที่ใช ้ในการ Service เป็น Random สมมุติว่าเป็น Exponential ด ้วยเวลาเฉลี่ย T

3.

ระบบสามารถรับ Customer ได ้ไม่จํากัด แต่จะ Service ได ้ครั้งละคน

4.

ถ ้าทุก Server เต็ม Customer ใหม่จะต ้องรอใน Queue ในกรณีนี้จะเกิด Queuing Delay 5.

ระบบนี้เรียก M/M/1 หรือ M/M/1/∞ แสดงได ้ด ้วย Simple Markov Model

6.

State Probability จะมีการกระจายแบบ Second Erlang (Erlang C) Distribution

0

1

2

i

j

p P = p P

i

ij

j

ji

Queuing System: Model

S

S

Service Rate

Arrival Rate = λ

= µ=1/Ts

Customer

S

Queue(FIFO)

Customer

S

M/M/N

Queuing System

Servers

M/M/1: Summary

λ

µ=1/Ts

S

Arrival = Poisson, λ

Inter Arrival Time = Exponential, 1/λ

Service Rate, µ

Service Time, Ts (1/µ) = Exponential

Queue = FIFO

1 Server

Queuing Model(1 Server);

M/M/1

Queue = 0, No Delay

Queue = Delay

0

1

N+1

N+2

X

Server ว่าง Server Busy

1/Ts = service rate

arrival rate

For each server

λ

µ =1/ Ts

การทํางานของ M/M/1

State = 0

No Q Delay

Queue Empty

State = 1

Delay

State = x; Queue = infinity

Customer Wait in Q

State = x; x = Q+1

Severe Delay

Queue Overflow (Full)

กรณีที่

Congestion

Queue มีขนาดจํากัด = Q

Packet Lost

เปรียบเทียบ Queuing Model

(N Server); M/M/N

Queue = 0, No Delay

Queue = Delay

0

1

i

j

N

N+1

N+2

X

Server ว่าง

i Server Busy

N Server Busy

1 Server Busy

1/h = service rate

A/h=arrival rate

For each server

Maximum Service Rate

= N/h

Service Rate at State

k = k/h

M/M/N

State < N

No Delay

Queue Empty

State >= N

Delay

Customer Wait in Q

Limited Q

Severe Delay

Queue Overflow (Full)

Congestion

Network Model using

M/M/1

แต่ละ Router เชื่อมต่อกันด ้วย Logical Link เดียว

เสมือนว่ามี Transmitter ตัวเดียวในการส่งข ้อมูลผ่าน Link

Network Model (M/M/1)

สมมุติว่าคอมพิวเตอร์ที่ต่อกับ Router

ต ้องการส่งข ้อมูลถึงกัน ตามเส ้นทางที่กําหนด

Network Model (M/M/1)

ที่ Output ของ Router สามารถ Model โดยใช ้ M/M/1

Network Model (M/M/1)

ถ ้าเราให ้ทุก Model เป็น M/M/1

ดังนั้น Delay จะเป็นผลรวมของ Delay แต่ละอัน

Kendal Notation

Kendal Notation

Kendal Notation

Next Week

• M/M/1 Analysis and Examples

• HW 5 Due

CPE 332

Computer Engineering

Mathematics II

Week 7

Part II, Chapter 6

Queuing System Continue

Extra: PRNG

Topics

• Single Server, M/M/1

• Kendal Notation

• Applications

Queuing System Case 3: Delay System

Limited Server=N; With Unlimited Queue

x

ρ

λ < µ

p ; 0 ≤ x N

1

0

N

1



k

!

x

Nρ

ρ 

k

−λ

λ

p =

p

x

x

=

+

N

N

e

; 0

ρ

[

P k] =

N

N!( N − ρ)

k

k =

!

t / T

0

s

1

  p ; x N

[

P T

t] = e

; T =

k!

s

0

µ

 N!  N

Arrival Rate = λ

Departure Rate = µ

Customer

System, N Server

Customer

1.

สมมุติว่า Customer แต่ละคนที่เข ้ามาเป็น Poisson และได ้รับการ Service จากระบบทันที

2.

เวลาที่ใช ้ในการ Service เป็น Random สมมุติว่าเป็น Exponential ด ้วยเวลาเฉล่ีย T

3.

ระบบสามารถรับ Customer ได ้ไม่จํากัด แต่จะ Service ได ้สูงสุด N พร ้อมๆกัน

4.

ถ ้าทุก Server เต็ม Customer ใหม่จะต ้องรอใน Queue ในกรณีนี้จะเกิด Queuing Delay 5.

ระบบนี้เรียก M/M/N หรือ M/M/N/∞ แสดงได ้ด ้วย Simple Markov Model

6.

State Probability จะมีการกระจายแบบ Second Erlang (Erlang C) Distribution

0

1

2

N

N+1

p P = p P

i

ij

j

ji

Queuing System Case 3: Delay System

Server=1; With Unlimited Queue; M/M/1

λ < µ 1

 ρ

−

p =

+1

=1− ρ

0

k

−λ

λ

 − ρ

e

1

(

)

[

P k] =

t / Ts

1

x

[

P T

t] = e

; T =

k!

s

p = ρ 1

( − ρ)

µ

x

Arrival Rate = λ

Departure Rate = µ

Customer

System, 1 Server

Customer

1.

สมมุติว่า Customer แต่ละคนที่เข ้ามาเป็น Poisson และได ้รับการ Service จากระบบทันที

2.

เวลาที่ใช ้ในการ Service เป็น Random สมมุติว่าเป็น Exponential ด ้วยเวลาเฉล่ีย T

3.

ระบบสามารถรับ Customer ได ้ไม่จํากัด แต่จะ Service ได ้ครั้งละคน

4.

ถ ้าทุก Server เต็ม Customer ใหม่จะต ้องรอใน Queue ในกรณีนี้จะเกิด Queuing Delay 5.

ระบบนี้เรียก M/M/1 หรือ M/M/1/∞ แสดงได ้ด ้วย Simple Markov Model

6.

State Probability จะมีการกระจายแบบ Second Erlang (Erlang C) Distribution

0

1

2

i

j

p P = p P

i

ij

j

ji

M/M/1: Summary

λ

µ=1/Ts

S

Arrival = Poisson, λ

Inter Arrival Time = Exponential, 1/λ

Service Rate, µ

Service Time, Ts (1/µ) = Exponential

Queue = FIFO

1 Server

Queuing Model(1 Server);

M/M/1

Queue = 0, No Delay

Queue = Delay

0

1

N+1

N+2

X

Server ว่าง Server Busy

1/Ts = service rate

arrival rate

For each server

λ

µ =1/ Ts

การทํางานของ M/M/1

State = 0

No Q Delay

Queue Empty

State = 1

Delay

State = x; Queue = infinity

Customer Wait in Q

State = x; x = Q+1

Severe Delay

Queue Overflow (Full)

กรณีที่

Congestion

Queue มีขนาดจํากัด = Q

Packet Lost

Kendal Notation

Kendal Notation

Kendal Notation

Analysis ય઼�શ઼ M/M/1

• สมมุติตอนแรกว่า Queue มีขนาดไม่จํากัด

• ใช ้ M/M/1 ในการ Model แต่ละ Port ของ Router (หรือ

Switch L3)

• Arrival คือจํานวน Packet ที่เข ้ามาในช่วงเวลาหนึ่ง ปกติ

วัดเป็น pps

• ขนาดของ Packet สมมุติว่าไม่แน่นอน แต่มีการกระจาย

แบบ Exponential

– Service Time ของแต่ละ Packet จะเป็น Exponential ด ้วย ทั้งนี้

ขึ้นอยู่กับค่า Link Speed ของ Output Port

• ค่า Server Utilization เท่ากับอัตราส่วนของ Arrival

Rate หารด ้วย Service Rate จะบ่งบอกอัตราส่วนที่

Server จะ Busy และคือ Link Utilization ของ Output

Port ด ้วย

ρ = λ / µ

Queuing in Communication

NW and M/M/1

Arrival Rate

Service Rate = 1/Service Time

λ

µ = 1/ Ts

ρ = λ / µ = λ Ts

Example

• Router ได ้รับ Packet เฉลี่ย 8 pps

– ความยาวของ Packet มีการกระจายแบบ Exponential

ด ้วยความยาวเฉลี่ย 500 Octet

– Link ที่จะส่งออกไป มีความเร็ว 64 kbps

• 1. Arrival Rate,λ = 8 pps

• 2. ความยาวเฉลี่ยของ Packet = 4000 bit

• 3. ความเร็ว Link = 64 k ดังนั้น Service Time, Ts

ของแต่ละ Packet = 4000/64k = 1/16

• 4. Service Rate(µ) = 16 pps

• 5. Server Utilization = 8/16 =0.5 = 50%

Assumption

• 1. อย่าลืมว่า Packet ที่เข ้ามา ต ้องเป็น

Independent และ Random มันจึงเป็น

Poisson

• 2. ความยาวของ Packet จะสมมุติว่าเป็น

Exponential ดังนั้น Service Time จะเป็น

Exponential ด ้วย แม ้ว่าสมมุติฐานนี้จะไม่

ถูกต ้องนัก

• 3. มี Output Link เดียว คือเป็น Single

Server

• 4. ดังนี้แล ้ว จึงจะเป็น M/M/1

Utilization

• Utilization บอกอัตราส่วนที่ Server จะ

Busy และสัมพันธ์กับ Probability ที่

Queue จะว่าง

– Probability ที่ Server ว่าง 1 − ρ

• ใน Network คือคือ Probability ที่

Output Link จะ Busy ด ้วย

ρ = λ / µ = λ Ts

Arrival Rate

• เนื่องจาก Arrival Rate มีการกระจายแบบ

Poisson ดังนั้นถ ้าให ้

λ เป็นอัตราเฉลี่ย

ของ Customer (Packet) ที่เข ้ามาใน

ช่วงเวลา 1 วินาที

– Probability ที่จะมี k customer (Packet) เข ้า

มาในช่วงเวลา T วินาทีจะหาได ้จาก

( T ) k −λ

λ

e T

p( X = k) = p( k) =

k!

−λ T

p( )

0 = e

Service Time

Ts เป็น Service Time เฉลี่ย และ Service

Rate หาได ้จาก µ =1/ Ts

• เนื่องจาก Service Time เป็น Random

Variable ที่มีการกระจายแบบ Exponential

ดังนั้น Probability ที่ Service Time จะมี

ค่าน ้อยกว่า T จะเป็น

T / Ts

p t

( < T ) = 1 − e

Queue Distribution

• การกระจายของ Customer (State Probability)

สามารถคํานวณได ้จาก Probability ที่ ระบบ จะมี

k Packet อยู่ดังนี้

k

p = p T )

k

0

s

• โดยที่ p0 คือ Probability ที่ ระบบ จะว่าง

p = 1− ρ = 1− λ T

0

s

• ดังนั้น เราได ้

k

k

p = 1

( − λ T )(λ T ) = 1

( − ρ )ρ

k

s

s

• กล่าวคือ การกระจายของ Customer ในระบบ

หรือค่า State Probability จะเป็น Geometric

Distribution

લ઼�

�હ્યસ઼દ્ઘ� Customer હ્વ�����, N

લ઼�

�લ઼�

�હ્યસ઼દ્ઘ�ય઼�શ઼ State, 𝐸𝐸[𝑋𝑋]

จาก 𝑃𝑃 𝑋𝑋 = 𝑥𝑥 = 𝑝𝑝𝑥𝑥 = 1 − 𝜌𝜌 𝜌𝜌𝑥𝑥 เราได ้

𝐸𝐸 𝑋𝑋 = ∑∞

𝑥𝑥=0 𝑥𝑥 𝑝𝑝𝑥𝑥 = (1 − 𝜌𝜌) ∑𝑥𝑥=0 𝑥𝑥𝜌𝜌𝑥𝑥

แต่จาก (CPE231) ∑∞𝑥𝑥=0 𝑥𝑥𝜌𝜌𝑥𝑥−1 = 1

(1−𝜌𝜌)2

ดังนั้น

𝑋𝑋 = 1 − 𝜌𝜌 𝜌𝜌 ∑∞𝑥𝑥=0 𝑥𝑥𝜌𝜌𝑥𝑥−1 = 𝜌𝜌

1−𝜌𝜌

Queuing Delay

• จาก Geometric Distribution ค่าเฉลี่ย

คือจํานวน Customer เฉลี่ย คือจํานวน

Packet เฉลี่ยในระบบ จะหาได ้จาก

ρ

N = 1− ρ

• ถ ้าคิดเฉพาะจํานวน Customer เฉลี่ยใน

Queue เราจะได ้ ρ

ρ2

N = N − ρ =

− ρ =

Q

1 − ρ

1 − ρ

Queuing Delay

• สําหรับ Network ค่าเฉลี่ย จํานวน Packet

เฉลี่ยใน ระบบ และใน Queue จะหาได ้จาก

ρ

ρ 2

N =

N =

1 − ρ

Q

1− ρ

• ถ ้าแต่ละ Packet ต ้องใช ้เวลาเฉลี่ยในการ

Service

T

ดังนั้น ค่า

s

Queuing Delay

จะเป็น

ρ T

ρ

W =

s

=

1 − ρ

µ − λ

System Delay

• แต่ละ Packet ต ้องใช ้เวลาเฉลี่ยในการ

Service

T

ดังนั้น ค่า

s

Queuing Delay

จะเป็น

ρ T

ρ

W =

s

=

1 − ρ

µ − λ

• และเวลาเฉลี่ยทั้งหมดที่ลูกค ้าจะต ้องรอใน

ระบบทั้งหมดจะเป็น

ρ

T = W + T =

+ 1 = 1

s

µ − λ µ µ − λ

Little’s Theorem

• ถ ้า T เป็นเวลาเฉลี่ยที่ลูกค ้าอยู่ในระบบ และ

λ เป็น Arrival Rate ดังนั้นจํานวนลูกค ้า

เฉลี่ยในระบบจะเท่ากับ

N = λ T

N = λ W

Q

สรุป M/M/1

สรุป M/M/1

สรุป M/M/1

M/M/1 Example 1

6 5 4 3 2 1

λ = ?

µ = ?

S

M/M/1 Example 1

M/M/1 Example 1

M/M/1 Example 1

M/M/1 Example 1

M/M/1 Example 2

λ = ?

µ = ?

M/M/1 Example 2

4.ธนาคารต้องจัดท่ีน่ังก่ีท่ี(รวมท่ีน่ังตอนใช้บริการ)เพ่ือจะ

ให้แน่ใจว่าอย่างน้อย 80% ของผู้ท่ีเข้ามาจะได้น่ัง

4.ธนาคารต้องจัดท่ีน่ังก่ีท่ี(รวมท่ีน่ังตอนใช้บริการ)เพ่ือจะ

ให้แน่ใจว่าอย่างน้อย 80% ของผู้ท่ีเข้ามาจะได้น่ัง

𝑥𝑥 = 9

M/M/1 Example 3

M/M/1 Example 3

M/M/1 Example 3

End of Chapter 6

• HW6 Due Monday Noon

– ส่งที่พี่หนึ่งเท่านั้น ก่อนเที่ยง จันทร์ 29 ก.พ.

– ใส่ตะกร ้า หน ้าโต๊ะ Counter อย่าส่งผิดที่

– เฉลยจะประกาศวันถัดไปบน Web

• Next

– Random Number Generator

– Preparation for Midterm Exam

Topics

• Random Number and Properties

• Random Number Test

• Pseudo Random Number Generator

MT Exam Preparation

Random Number

• Random Number เป็น Set ของตัวเลขที่

สุ่มขึ้นมาแบบอิสระ(Random) และไม่ขึ้นต่อ

กัน (Independent) โดยมีค่าการกระจาย

(Distribution) ของชุดตัวเลข ตามที่

กําหนดแบบใดแบบหนึ่ง

• เลขแต่ละตัวที่สุ่มจะต ้องไม่มีความสัมพันธ์

กัน ทดสอบได ้จากค่า Correlation

મ઼��હ્વહ઼�

શ઼�� Random Number

• Gambling

• Statistical Sampling

• Simulation

• Cryptography

• Computer Games

• Hash Algorithm/Searching Algorithm

Random Number Generation

• เป็นอุปกรณ์ที่ใช ้เป็นตัวกําเนิดชุดของตัวเลข

Random ที่ต ้องการ

– True Random Number Generator(TRNG)

• ประกอบด ้วยอุปกรณ์ทาง Physic ที่ใช ้กําเนิดตัวเลข

• มักจะได ้จากแหล่งกําเนิดของ Noise ที่เกิดตาม

ธรรมชาติ

– Thermal Noise/Shot Noise/Radioactive

• หรือใช ้ Roulette Wheel

– Pseudo Random Number Generator(PRNG)

• Computational Algorithm จากสมการทาง

คณิตศาสตร์ โดยเริ่มจากตัวเลข “Seed”

• มีคุณสมบติเหมือนกับเป็น Random แต่มี Period และ

สามารถคํานวณล่วงหน ้าตามสมการได ้

PRNG: Pseudo Random Number Generator:

Linear Congruential Generator

• เริ่มจาก “Seed”, 𝑋𝑋0 และคํานวณ Series

𝑋𝑋1, 𝑋𝑋2, 𝑋𝑋3, … จากสมการ Recurrence

𝑋𝑋𝑛𝑛+1 = 𝑎𝑎𝑋𝑋𝑛𝑛 + 𝑏𝑏 𝑚𝑚𝑚𝑚𝑚𝑚 𝑚𝑚

– a, b และ m เป็น Integer ปกติจะมีขนาดใหญ่; จํานวนเลขสูงสุด

ที่ไม่ซํ้ากันจะถูกกําหนดด ้วยค่า Modulus m

• Sequence ที่ได ้จะมีค่าระหว่าง [0,m) หรือ 0 ≤ 𝑋𝑋 < 𝑚𝑚

• Seed 𝑋𝑋0; 0 ≤ 𝑋𝑋0 < 𝑚𝑚

• Multiplier 𝑎𝑎; 0 < 𝑎𝑎 < 𝑚𝑚

• Increment 𝑏𝑏; 0 ≤ 𝑏𝑏 < 𝑚𝑚

• ถ ้า b = 0 เราเรียก Multiplicative Congruential Generator หรือ Lehmer

Generator มิฉะนั้นแล ้วเราเรียก Mixed Congruential Generator

– ตัวเลขที่ได ้จะมีการกระจายแบบ Uniform ถ ้าต ้องการการกระจาย

แบบอื่น จะใช ้การ Transformation

– Algorithm นี้เรียก LCG หรือ Linear Congruential

Generator เป็นวิธีที่ง่ายในการสร ้าง PRN

Period of LCG

• Period มีค่าได ้สูงสุดคือ m เรียก Full

Period

– บางกรณีจะได ้น้อยกว่านั้น ขึ้นอยู่กับการเลือกค่า

‘a’ และในกรณีที่ b=0

• Generator จะมี Full Period ก็ต่อเมื่อ

(Hull-Dobell Theorem)

– 1. b และ m เป็น Relative Prime

– 2. a-1 จะต ้องสามารถหารได ้ด ้วยทุกๆ factor

ที่เป็นค่า prime ของ m

– 3. a-1 จะต ้องหารได ้ด ้วย 4 ถ ้า m หารได ้ด ้วย 4

Parameters used

• ส่วนใหญ่จะใช ้ LCG ที่ m มีค่าเป็นกําลังของ

2 ที่นิยมมากสุดคือ 232 และ 264

– MS Visual C++ ใช ้ m=232,

a=214013,b=2531011

– MS Visual Basic 6 ใช ้ m=224,

a=1140671485,b=12820163

લ઼�

ઙ્ખ દ્ભઢ્ઢ��

ઙ્ઘ�

LCG

• มีข ้อดีคือคํานวณได ้ง่าย แต่ให ้ Series ของ

Random Number ที่มีคุณสมบัติพอประมาณ

เท่านั้น เพราะให ้ค่า Serial Correlation สูง

– ไม่เหมาะกับการนําไปใช ้ใน Monte Carlo

Simulation

– ไม่เหมาะกับการนําไปใช ้ใน Cryptography

• สามารถสร ้างได ้จากวงจร Linear Feedback

Shift Register(LFSR)

– ปกติวงจรนี้จะผลิต Stream ของ bit

• ได ้ Uniform Distributed PRNG

PRNG ��� હ્યહ઼�

• Blum Blum Shub

• Wichmann-Hill

• Multiply with Carry

• Mersenne Twister (นิยมมากสุด)

• Park-Miller

• RC4

• Rule 30

Randomness Test

• เพื่อหา Pattern ในชุด หรือ Series ของตัวเลข

– ซึ่งไม่ควรจะมีอยู่ถ ้าชุดตัวเลขเป็น Random และ

Independent อย่างแท ้จริง

• ใช ้ในการทดสอบ และเปรียบเทียบ Algorithm ต่างๆ

ที่ใช ้ผลิต PRNG โดยใช ้พื้นฐานจาก

– Statistical Test ทดสอบจากคุณสมบัติทางสถิติ

Diehard Tests ประกอบด ้วยชุดของ Test

TestU01 library ประกอบด ้วย utilities สําหรับ

statistical testing ของ uniform RNG

– Transform ดูคุณสมบัติเมื่อ Transform

• Hadamard Transform(Generalized FT)

– Complexity ความซับซ ้อนของตัวเลข

• Kolmogorov complexity

Midterm Preparation

• เก็บ 35 เปอร์เซ็นต์

• ไม่มีการ Make Up

• สอบ 3 ชั่วโมง บทที่ 1-6

– Monday 7 March; 13:30 – 16:30

• ต ้องใช ้เครื่องคิดเลข

– รุ่นตามที่คณะประกาศเท่านั้น ห ้ามใช ้อุปกรณ์อื่นๆ

• สมการที่สําคัญจะให ้มา ไม่ต ้องจํา

6 ข้อ บทละ 1 ข้อ รวม 6 ข้อ เลือกทํา 5 ข้อ

5 ข้อ 50 คะแนน คิดเป็น 35 %

สมการที่ให้

ในการสอบ MT

Review દ્મ�

દ્ધય઼�

���મ઼દ્ભ��

• Chapter 1: Vector

– Magnitude/Direction/Unit Vector

– Direction Cosine

– Component Vector/Position Vector

– Dot/Cross Product and Properties

– Equation of Line and Plane, Angle of Vectors

• Chapter 2: Matrices

– Types of Matrices, Minor, Cofactor, Diagonal

– Determinant by Expansion

– Inverse of Matrix

– Rank/Reduced Matrix (Process of Elimination)

– Homogeneous/Non Homogeneous Linear Eq.

Review દ્મ�

દ્ધય઼�

���મ઼દ્ભ��

• Chapter 3: Eigen Value/Vector/Diag

– Eigenvalues

– Eigenvectors

– Diagonalization

– Symmetric/Orthogonal Matrix

• Chapter 4: Probability

– Conditional Probability/Bayes Rule

– Random Variable

– CDF/PDF/PMF; Poisson/Exponential Distribution

– Expectation Concept

– Mean, Mean Square, Variance

– Joint Random Variable, Correlation/Covariance

Review દ્મ�

દ્ધય઼�

���મ઼દ્ભ��

• Chapter 5: Random Process

– Stationary/Ergodic

– Autocorrelation/Cross Correlation

– Counting/Poisson Process/Birth and Death Process

– MarKov Model; Global Balanced Equation

– Simple MarKov Model; Detail Balanced Equation

• Chapter 6: Queuing System

– M/M/1 Concept

– Arrival Rate/Inter Arrival Time

– Departure Rate/Service Time

– Utilization, P[X=k], P[X<=k]

– Average Customer(System/Q), Time(System/Q)

CPE 332

Computer Engineering

Mathematics II

Week 10

Part III, Chapter 8

Error, Transcendental Function and Power Series

Taylor Series and Function Approximation

Today Topics

• Numerical Methods

• Errors

• Transcendental Function

• Series

• Power Series

• Taylor Series

• MacLaurin Series

• Function Approximation

• HW VII

• MT Exam and Cumulative Grades

Numerical Methods

• ในวิชาคณิตศาสตร์ที่เราเรียนที่ผ่านมา การหา

Solution จะประกอบด ้วยการแก ้สมการ ซึ่งรวมถึง

การแก ้สมการพิคณิต การหา Derivative การ

Integrate การแก ้สมการ Differential

Equation รวมถึงการใช ้เครื่องมือทางคณิตศาสตร์

เพื่อแก ้สมการ เช่น Fourier Transform,

Laplace Transform หรือ Z-Transform

วิธีการดังกล่าวจัดว่าเป็นการแก ้ปัญหาที่เรียก

Analytical Method อย่างไรก็ตาม ไม่ใช่ว่า

สมการทุกสมการจะหาคําตอบได ้(ที่อยู่ในรูปของ

Explicit Form) ดังนั้นเราจําเป็นที่จะต ้องหาวิธีอื่น

ในการแก ้ปัญหา นั่นก็คือวิธีการที่เรียก Numerical

Method

Numerical Methods

• วิธีการทาง Numerical ที่เคยใช ้กันมาได ้แก่วิธีการ

Graphical Method และใช ้บวนการ

Interpolation แต่เมื่อมีการประดิษฐ์คอมพิวเตอร์

ขึ้นมา ก็มีการนําคอมพิวเตอร์เข ้ามาใช ้ และขยาย

ขอบเขตของวิชานี้ออกเป็นสาขาใหม่ทาง

คณิตศาสตร์ มีการศึกษา Algorithm และวิธีการ

เขียนโปรแกรม รวมถึงการวิเคราะห์ค่าความ

ผิดพลาด(Error Analysis) และลักษณะการ

Convergence และ Complexity ของ

Algorithm

Y=f(x)=2x3-4x2-7x+5

• หา f(k) เช่น f(0), f(-2), f(7) ง่าย

• แก ้สมการ f(x)=k จะยาก เช่น หาค่า x ที่ทําให ้ f(x) = 0

Quadratic Equation

• Polynomial Degree = 1, Straight Line

– y = f(x) = mx+b; x = (y-b)/m

• Polynomial Degree = 2

• y = f(x) = ax2+bx+c

• แก ้สมการ y = f(x) = 3 เราได ้

– y’ = f’(x) = ax2+bx+(c-3) = 0

– x = [-b±√(b2-4a(c-3))]/2a

• Polynimial Degree = n, n > 2

– ไม่มีสมการโดยตรงในการแก ้,

– แต่มี Algorithm ในวิธีของ Numerical Method

• Algorithm จะเป็น Iterative Method

– คําตอบของ Numerical Method จะให ้แค่ค่าประมาณ

• ค่าจะถูกต ้องขึ้นเรื่อยๆ ใน Iteration ที่สูงขึ้น (คํานวณนานขึ้น) ถ ้า

Algorithm Converge

Start

Iterative

Iteration = 0

Set et

Technique

Init. X0

[Iteration i]

Iteration += 1

โปรแกรมจะ Converge

Calculate xi

ถ้า ei < ei-1

Estimate ei

N

Y

Error e

Stop

i < es

Errors

• ในการนําคอมพิวเตอร์มาช่วยแก ้ปัญหา

ทางด ้านคณิตศาสตร์นั้น จะต ้องเข ้าใจก่อน

ว่าลักษณะการทํางานของคอมพิวเตอร์ หรือ

การคํานวณจะขึ้นอยู่กับวิธีการที่เราเขียน

โปรแกรมและ Algorithm ที่ใช ้ ในการ

คํานวณโดยใช ้เครื่องคํานวณใดใดจะมี

Error เกิดขึ้นเสมอ ซึ่งจะแบ่งได ้เป็นสอง

ประเภทคือ Round- Off Error และ

Truncation Error

Errors

• Round- Off Error เก ิดจากการเก็บ

ต ัวเลขใน Memory ของคอมพ ิวเตอร ์น ั�น

จะใช ้ขนาดของ Memory ที่จําก ัด และ

ข �ึนอยู่ก ับว่าเราให้ต ัวแปรชน ิดอะไร เช ่น

Integer, Long Integer, Float, หร ือ

Double

Errors

• Truncation Error เก ิดจาก Algorithm

หร ือวิธีที่เราให้คอมพ ิวเตอร ์คํานวณ น ั�น

เป็นแค่การประมาณ จากสมการทาง

คณ ิตศาสตร ์ที่เราต ้องการหาค่าโดยใช ้

คอมพ ิวเตอร ์เข ้าช ่วย และเราไม่ได ้แก ้

สมการทางคณ ิตศาสตร ์จร ิงๆ

– การคํานวณ จะต ้องม ีการบวกก ันของ

Infinite Terms จึงจะได ้ค่าที่ถูกต ้อง

ด ังน ั�นคอมพ ิวเตอร ์จะคํานวณ ค่าประมาณ

โดยต ัดเทอมท้ายๆทิ�ง

Precision, Accuracy,

Significant Digit

• ดังนั้นเราจําเป็นที่จะต ้องรู ้ว่า Error เรามีมากแค่ไหน และ

พอจะยอมรับได ้หรือไม่ นั่นขึ้นอยู่กับการนําการคํานวณไป

ใช ้งาน ปกติตัวเลขที่จะไปใช ้งานทางวิศวกรรมศาสตร์จะ

กําหนดด ้วยค่า Significant Digit และจะยอมให ้ความไม่

แน่นอนของตกเลขตกอยู่ที่หลักท ้ายเท่านั้น ดังนั้นค่าของ

Significant Digit จะบ่งบอกถึงจํานวนหลักของตัวเลขที่

จะนําไปใช ้ได ้อย่างมั่นใจ

• คําว่า Accuracy หมายถึงตัวเลขที่เราได ้นั้นมีค่าใกล ้เคียง

กับค่าในความเป็นจริงเท่าไร แต่คําว่า Precision หมายถึง

ตัวเลขที่ได ้มานั้นเกาะกลุ่มกันมากแค่ไหน ดังนั้นค่าของ

Precision จะเป็นตัวกําหนดจํานวนของ Significant

Digit ที่จะต ้องใช ้ และค่าของ Error จะเป็นตัวกําหนด

ความ Accuracy ของตัวเลข

การยิงเป้า

Precision = ตํ่า

Precision = สูง

Accuracy = ตํ่า

Accuracy = ตํ่า

Precision = ตํ่า

Precision = สูง

Accuracy = สูง

Accuracy = สูง

Accuracy ว ัดจากค่า Error ที่เกิดระหว่างค่าเฉล่ียของ Data ก ับค่าจริง

Precision ว ัดจากค่า Variance หรือ SD ของกลุ่มของ Data

Significant Digits (Figures)

• คือตัวเลขที่มีความหมายในการกําหนดค่า

Precision

– เป็นเลขทุกตัวยกเว ้นเลขศูนย์นําหน ้า และศูนย์ที่

ต่อท ้าย

• จํานวนศูนย์ต่อท ้ายที่เป็นส่วนของ s.f. บางทีเรากําหนด

โดยใช ้เครื่องหมาย ‘bar’

• Examples

– 11.152 = 5 s.f.

– 0.00879 = 3 s.f.

– 000125600. = 6 s.f.

– 0123400000 = 7 s.f.

• เลขทางขวามือของ s.f. จะไม่นํามาใส่ ให ้ทํา

การ Round-Off เช่น

– 456.389864 (6 s.f.) = 456.390

– 1287563.94 (3 s.f.) = 1290000

Error

• คือค่าที่แตกต่างจากค่าที่แท ้จริง

• เป็นตัวกําหนดค่า Accuracy

• แบ่งเป็น

– Absolute Error, Et

– Relative Error (%จากค่าจริง), et

• นอกจากนี้ยังมี

– Estimated Error, ea

• ค่าประมาณของ Error (Relative Error)

Absolute Error

Relative Error

Error Estimation

• Convergence

• Divergence

Mean Square Error(MSE)

• ในการเปรียบเทียบความแตกต่างระหว่างค่าที่

แท ้จริงกับค่าที่ประมาณได ้ สําหรับกลุ่มของ

ตัวอย่าง เรามักจะใช ้ค่าเฉลี่ยของ Error ในรูปของ

ค่าเฉลี่ยของกําลังสองของ Error ในแต่ละคู่ เรียก

Mean Square Error

– ถ ้าให ้ Y เป็นค่าที่แท ้จริง และ

Yˆ เป็นค่าที่ประมาณได ้

สําหรับคู่ของตัวอย่าง N ตัวอย่าง

1

1 N

MSE =

∑− ˆ

( Y

2

Y )

i

i

N i=0

– ค่า RMSE หรือค่า Root Mean Square Error คือค่า

Square Root ของค่า MSE

N 1

RMSE =

1 ∑− ˆ( Y − 2

Y )

i

i

N i=0

Functions

• เป็นความสัมพันธ์แสดง Set ของ Input

และ Set ของ Output

– Function สามารถรับค่า Input ได ้หลาย Set

(หลายตัวแปร)

– การ Mapping เป็นทางเดียว

• ถ ้า Map กลับจะเป็น อีก Function ที่เป็น Inverse

Function กับตัวเดิม

• Inverse อาจจะไม่มีสําหรับบาง Function

x

y

f(.)

y=f(x)

x=f-1(y) อาจจะไม่มี

Polynomial Functions

• เป็นรูปแบบของ Function ที่ง่ายที่สุด

• Polynomial function Degree ‘n’ จะอยู่ใน

รูปของ(One Argument)

n

n

y = f ( x) = a x + a

x 1 +  + a x 2 + a x + a

n

n −1

2

1

0

n

n

= ∑ i

a x = a +

a x

i

0

ii

i =0

i =1

Polynomial of degree 2:

f(x) = x2 - x - 2

= (x+1)(x-2) Polynomial of degree 3:

f(x) = x3/4 + 3x2/4 - 3x/2 - 2

= 1/4 (x+4)(x+1)(x-2)

Polynomial of degree 4:

f(x) = 1/14 (x+4)(x+1)(x-1)(x-3) + 0.5

Polynomial of degree 5:

f(x) = 1/20 (x+4)(x+2)(x+1)(x-1)(x-3) + 2

เอามาจาก

http://en.wikipedia.org/wiki/Polynomial

Polynomial of degree 6:

Polynomial of degree 7:

f(x) = 1/30 (x+3.5)(x+2)(x+1)(x-1)(x-3)(x-4) + 2

f(x) = (x-3)(x-2)(x-1)(x)(x+1)(x+2)(x+3)

Transcendental Function

• เป็น Function ที่ไม่ใช่ Algebraic

– Algebraic function คือ Function ที่สามารถ

นิยามได ้จาก Root ของสมการ Polynomial

• Transcendental Function ไม่สามารถ

เขียนในรูป Solution ของ Polynomial

– Exponential Function

– Logarithm

– Trigonometric Functions

• ด ้วยเหตุผลนี้ การหาค่าของ Function

อาจจะต ้องใช ้วิธีการประมาณค่าแทน

Transcendental Examples

f ( x) = π

x

1

f ( x) = cx , c ≠ 1

,

0

2

f ( x) = xx

3

1

f ( x) = x x

4

f ( x) = log x, c ≠ 1

,

0

5

c

f ( x) = sin x

6

การประมาณค่าของ Function

• ในส่วนนี้เราจะ Review เรื่อง

– Infinite Series

– Power Series

– Taylor Series

– Maclaurin Series

– การประมาณค่าโดยใช ้ Series

Series และ Infinite Series

• Series คือผลรวมของเทอมใน Sequence ใน

กรณีของ Infinite Series นั้น เทอมที่จะต ้อง

นํามาบวกกันจะไม่มีที่สิ้นสุด

• ถ ้าให ้แต่ละเทอมใน Sequence เป็น an ดังนั้น

Series คือ a1+a2+a3+… เช่น

S = ∑ a = a + a + a +

n

1

2

3

n=1

a = 1

S =

n

n

∑ 1

;

= 1 + 1 + 1 +

2

n

n=1 2

2

4

8

• NOTE: บางครั้ง Index ของ Summation อาจจะเริ่มที่ 0 ก็ได ้

Convergence of Series

• Series จะ Converge ก็ต่อเมื่อผลรวมของ

Term สามารถหาค่าได ้

– กล่าวคือผลรวมของเทอมจนถึง N มี Limit

N

S = ∑

a = lim S = lim

a

n

N

N →∞

N →∞

n

n=0

n=0

• ถ ้า Series ไม่ Converge มันจะ Diverge

Power Series

• Power Series เป็น Infinite Series ใน

รูปของ (c = Constant, an = coefficient,

x เป็น Variable ที่อยู่รอบๆ c)

f ( x) = ∑

a ( x

n

c)

n

n=0

– Power Series ที่สําคัญคือ Taylor Series

• ในกรณีที่ c=0 เราได ้

f ( x) = ∑∞ a xn = a + a x+ 2

a x +

3

a x + 

n

0

1

2

3

n=0

– เช่นในกรณีของ Maclaurin Series

Taylor Series/Maclaurin Series

• ให ้ f(x) เป็น Function ที่สามารถหา

Derivative ได ้อย่างไม่จํากัดในช่วง

ใกล ้เคียงกับ a ดังนั้น Taylor Series ของ

f(x) คือ Power Series ในรูปของ

(3)

( n)

f a

f

a

f

a

f

a

f ( a) +

'( )

''( )

( )

( )

2

3

n

( xa)+

( x a) +

( x a) + = ∑

( x a)

!

1

!

2

!

3

n

n=

!

0

• ในกรณีที่ a = 0 บางครั้งเราเรียก Maclaurin

Series

(3)

( n)

f

f

f

f

f ( x) = f ( )

0 +

'( )

0 x + ''( )

0

(0)

(0)

2

x +

3

x + = ∑

n

x

!

1

!

2

!

3

n

n=

!

0

Example 1

(3)

(4)

f ''( a)

f

( a)

f

( a)

• จาก

f ( x) = f( a)+ f '( a)( xa)+

( x

2

a) +

( x

3

a) +

( x

4

a) + ...

!

2

!

3

!

4

(3)

(4)

f

f

f

f ( x) = f (0) + f '(0) x +

''( )

0

( )

0

( )

0

2

x +

3

x +

4

x + ..

!

2

!

3

!

4

Suppos

e f ( x) = x

e , f '( x) = x

( n)

e ,..., f

( x) = x

e

( n)

f

( )

0 = 1

2

3

4

k

x

x

x

x

ha

We

ve f ( x) = x

e = 1+ x +

+

+

+ ... = ∑∞

!

2

!

3

!

4

k

k =

!

0

Example f

ind f ( )

2

2

= e

2

4

Second O

rder A

pproximat ion : e ≈ 1 + 2 +

= 5

2

2

4

8

16

Forth Order A

pproximat ion : e ≈ 1 + 2 +

+ +

= 7

2

6

24

2

4

8

16

32

64

Sixth Order A

pproximat ion : e ≈ 1 + 2 +

+ +

+

+

= .

7 356

2

6

24

120

720

2

4

8

16

32

64

128

256

Eighth Order A

pproximat ion : e ≈ 1 + 2 +

+ +

+

+

+

+

= 3873

.

7

2

6

24

120

720

5040

40320

Actual V

alue : 2

e = 7 389056099

.

นักศึกษาหาค่าประมาณของ e-2 ถึง 4 Significant Digit ได ้หรือไม่

Example 2:

(3)

(4)

f

f

f

f ( x) = f ( )0 + f '( a) x + ''( )0

( )

0

( )

0

2

x +

3

x +

4

x + ..

!

2

!

3

!

4

Suppos

e f ( x) = sin x, f '( x) = cos x, f ''( x) = − sin x, f '''( x) = −

(4)

cos x, f

= sin( x),...

3

5

7

n

x

x

x

ha

We

ve f ( x) = sin x = 0 + x + 0 −

+ 0 +

+ 0 −

+ ... = ∑∞ (− )1

2 n+1

x

!

3

!

5

!

7

n

n=

(2

)!

1

0

+

Example f

ind f (0 )

5

.

= sin .

0 5 radian

Second O

rder A

pproximat ion s

: in 0 5

.

5

.

0

5

.

0 3

Fort O

h rder A

pproximat ion : sin 0 5

.

5

.

0

= 47916

.

0

667

!

3

5

.

0 3

5

.

0 5

Sixth Order A

pproximat ion : sin 0 5

.

5

.

0

+

= 47942

.

0

708

!

3

!

5

5

.

0 3

5

.

0 5

5

.

0 7

Eighth Order A

pproximat ion : sin 0 5

.

5

.

0

+

= 47942

.

0

5533

!

3

!

5

!

7

Actual V

alue : sin 0 5

. = .

0 479425538

นักศึกษาหาค่าประมาณของ cos 0.5 ถึง 8 Significant Digit ได ้หรือไม่

Taylor (Maclaurin) Series

n

2

3

4

x

x

x

x

xe =∑ =1+ x+ + + +;∀ x

n

n=

!

!

2

!

3

!

4

0

n

x

log(1− x) = −∑

; −1 ≤ x < 1

n

n=1

n

x

log(1+ x) = ∑

n+

(−

1

)

1

; −1 < x ≤ 1

n

n=1

1

= ∑∞ n

x ; x < 1

1− x

n=0

m+

1−

1

m

x

= ∑ n

x ; x ≠ 1

1− x

n=0

x

=

n

nx ; x

1

2

∑∞

<

1

( − x)

n=0

n

n

1+ x

∑∞ (−

=

)

1 (2 )!

n

x = 1+ 1 x − 1 2

x + 1 3

x − 5

4

x + ; x < 1

2

2

8

16

128

n n

n

n=

1

(

2 )( !) (4 )

0

n

3

5

x

x

sin x = ∑

(− )

1

2 n+1

x

= x

+

−; ∀ x

n

n=

(2

)

1 !

!

3

!

5

0

+

n

x

x

cos x

∑∞ (−

2

4

=

)

1

2 n

x

= 1−

+

−; ∀ x

n

n=

(2 )!

!

2

!

4

0

MATLAB Program

• จะสาธิตการเขียน Function โดยใช ้ MATLAB

– เขียน Function ที่รับค่า Vector x และ Vector y

return z1 และ z2; x=-3:.1:3; y = -3:.1:3

– คํานวณค่า z1 = 2sin(x^2+y^2)/(x^2+y^2)

– คํานวณค่า z2 = xsiny-ycosx

– mesh (x,y,z1), (x,y,z2), (x,z1), (y,z2) และ (z1,z2)

– Surf function, shading modes

– Contour plot

• ศึกษาวิธีการเขียน Function ใน MATLAB โดยดู

จาก Tutorial 4-5

END OF WEEK 10

Download HW 7

Next Week Chapter 9: Zeros of

Functions

MATLAB Program

function [z1,z2]=test(x,y)

% function [z1,z2]=test(x,y)

% Test matlab program calculate z1=f(x,y)=sin(x^2+y^2)/(x^2+y^2)

% and z2=f(x,y)=xsiny-ysinx

n=length(x);

m=length(y);

z1=zeros(n,m);

z2=zeros(n,m);

for i = 1:n

for j= 1:m

t=x(i)^2+y(j)^2;

if (t ~= 0)

z1(i,j)=sin(t)/t;

else

z1(i,j)=1;

end

z2(i,j)=x(i)*sin(y(j))-y(j)*cos(x(i));

end

end

MATLAB Program

• function view2(x,y,z,az)

• mesh(x,y,z);

• view(0,az);

• for i = 0:10:360

• view(i,az);

• pause(0.05);

• end

MATLAB Program

• function view3(az)

• view(0,az);

• for i = 0:10:360

• view(i,az);

• pause(0.05);

• end

CPE 332

Computer Engineering

Mathematics II

Week 11

Part III, Chapter 9

Roots of Equations

Today Topics

• Roots of Equation

• Bracketing Method

– Bisection Method

– False Position Method

• Open Method

– Simple One-Point Iteration

– Newton-Ralphson Method

– Secant Method

– Multiple Roots Problem

– Modified Method

Roots of Equations

Zeroes of Functions

• กําหนด Function y=f(x)

– Zeroes ของ Function คือค่าของ x ที่ทําให ้

f(x)=0 หรือค่า y=0

– คือราก(root) ของสมการ f(x)=0 นั่นเอง

• ตัวอย่าง y=x3-2x2+x-5

– สมการ Polynomial, degree = 3

– ถ ้าเกิน Quadratic (Degree = 2) การหาราก

จะทําได ้ยาก, ไม่มีวิธีทาง Analytic Method

โดยตรง เหมือน y=Ax2+Bx+C=0

2

b ± b − 4 ac

x =

2 a

Zeroes of Functions

• ตัวอย่าง y=x3-2x2+x-5

– วิธีทาง Numerical ที่เราเคยเรียนมาในวิชา

เรขาคณิตย์ คือ Plot Graph และหารากของ

สมการจาก Graph

• 1. แต่วิธีนี้จะให ้ค่าที่หยาบ

• 2. กรณีที่เป็น Complex Root จะหาไม่ได ้ จาก

Graph

Zeroes of Functions

y=x3-2x2+x-5

2.4

Zeroes of Functions

y=x3-2x2+x-5

4334

.

2

x = − 2167

.

+ j 417

.

1

− 2167

.

j 417

.

1

2.4

Roots of Equation

• วิธีที่จะกล่าวต่อไป

– Numerical Method = Algorithm

– หา Zero Crossing (รากที่เป็นค่าจริง)

– สามารถได ้คําตอบให ้ถูกต ้องด ้วย Significant

Digit มากเท่าที่เราต ้องการ

• ต ้อง Run Algorithm นานขึ้น

• Algorithm จะต ้อง Converge

– เป็น Iterative Method

• ถ ้า Algorithm Converge, แต่ละ Iteration ที่

Run จะให ้คําตอบที่ถูกต ้องขึ้นเรื่อยๆ

• จะช ้าหรือเร็วขึ้นอยู่กับ Convergence Rate ของแต่

ละ Algorithm

Methods

• Bracketing Method

– Bisection

– False Position

• Open Method

– Simple one-point Iteration

– Newton-Ralphson

– Secant Method

Bracketing Method

• ใช ้หลักความจริงที่ว่า เมื่อ f(x)=0 มันจะต ้อง

เปลี่ยนเครื่องหมาย

– ดังนั้น f(x-)f(x+) < 0 เมื่อ f(x)=0

• เราเริ่มจาก สองค่าของ x คือ xl และ xu ที่มี

คุณสมบัติ f(xl)f(xu) < 0

– อย่างน ้อยต ้องมีคําตอบหนึ่งอยู่ในช่วงนี้

• Algorithm จะ Search หาคําตอบในช่วง

Bracket นี้ โดยจะลดขนาดของ Bracket ลง

เรื่อยๆ

Bracketing Method

f ( x ) f ( x ) < 0

l

u

( x , f ( x ))

u

u

xl

xu

( x , f ( x ))

l

l

Bracketing Method

f ( x ) f ( x ) = 81

− < 0

l

u

(

)

9

,

7

x = 4

l

x = 7

u

( ,

4 − )

9

Bracketing Method:Bisection

Bracketing Method

x = (4 + 7) / 2 = 5.5

r

(

)

9

,

7

x = 4

l

x = 7

u

( ,

4 9

− )

Bracketing Method

x = (4 + 7) / 2 = 5.5

r

(

)

9

,

7

x = 4

x = 7

l

u

( x , f ( x ))

r

r

.

5

(

,

5 − .

2

)

25

( ,

4 9

− )

xu=xr

Bracketing Method

(

)

9

,

7

x = 7

u

x = x =

5

.

5

l

r

Bracketing Method

(

)

9

,

7

x = 7

u

x = x =

5

.

5

l

r

x =

5

.

5

(

+ 7) / 2 = 25

.

6

r

Bracketing Method

x = x =

5

.

5

x = x

l

r

u

r

x =

5

.

5

(

+ 7) / 2 = 25

.

6

r

Bracketing Method:Bisection

Bracketing Method:Bisection

Bracketing Method:Bisection

.

667 38

f ( c) =

1

(

10

c /

1

.

68

e

) − 40 = 0

c

f

)

12

(

= 0669

.

6

38

.

667

f ( c) =

1

(

10

c /

1

.

68

e

) − 40 = 0

c

xu

x

x

l

r

f

)

16

(

= 2

− 2688

.

x = 12

(

+ )

16 / 2 = 14

r

f x

=

f

)

14

(

= 5687

.

1

(

)

1 5687

.

r

(2)

f

)

12

(

= 0669

.

6

x

= x = 14

u

r

38

.

667

f ( c) =

1

(

10

c /

1

.

68

e

) − 40 = 0

c

xu

x

x

l

r

f

)

16

(

= 2

− 2688

.

x = 12

(

+ )

16 / 2 = 14

r

f

)

14

(

= 5687

.

1

x = 14

(

+ )

16 / 2 = 15

38

.

667

r

f ( c) =

1

(

10

c /

1

.

68

e

) − 40 = 0

c

x

x

r

u

xl

f

)

16

(

= 2

− 2688

.

Bracketing Method: Bisection

False-Position Method

False-Position Method

False-Position Method

False-Position Method

False-Position Method

False-Position Method

Open Method:

Simple One-Point Iteration

Open Method:

Simple One-Point Iteration

Open Method:

Simple One-Point Iteration

Open Method:

Simple One-Point Iteration

Open Method:

Simple One-Point Iteration

Open Method:

Newton-Ralphson Method

Open Method:

Newton-Ralphson Method

Open Method:

Newton-Ralphson Method

Open Method:

Newton-Ralphson Method

Open Method:

Newton-Ralphson Method

Open Method:

Secant Method

Open Method:

Secant Method

Open Method:

Secant Method

Multiple Roots

Multiple Roots

Multiple Roots

Multiple Roots

Multiple Roots

Multiple Roots

Multiple Roots

Comparison

Comparison

Comparison

Homework 8

• นักศึกษาต ้องเขียนโปรแกรมช่วยคํานวณ

หรือใช ้ Spreadsheet (MS Excel) ช่วย

คํานวณ

– แนะนําให ้ใช ้ MATLAB

• เขียน Function หรือ Scratch File ก็ได ้

• หรือคํานวณจาก Workspace โดยตรง

• Download คําถามและตอบคําถาม

CPE 332

Computer Engineering

Mathematics II

Week 12

Part III, Chapter 10

Linear Equations

Today(Week 13) Topics

• Chapter 9 Linear Equations

– Gauss Elimination

– Gauss-Jordan

– Gauss-Seidel

– LU Decomposition

– Crout Decomposition

• HW(9) Ch 9 Due Next Week

MATLAB Programming

• เราสามารถเขียน Function การคํานวณโดยใช ้ MATLAB

Editor และบันทึกเป็น ‘.m’ File

– ขึ้นบันทัดแรกของ Function ด ้วย

function [List ของค่าที่ส่งคืน]=fname(List ของ Parameter)

function [x,y,z]=find123(a,b,c)

– ภายใน Function สามารถใช ้ Loop, Branch ได ้เหมือนการ

เขียนโปรแกรม, สามารถกําหนด Local Variable ภายในได ้

เช่นกัน

– อย่าลืมว่า พื้นฐาน Variable จะเป็น Matrix

• Function นี้สามารถเรียกใช ้งานได ้ใน MATLAB

• ดูรายละเอียดใน Tutorial 4-5 ของ MATLAB

Ex: หารากของ

− cos x

− sin x

f ( x) = sin 3 x e

+ cos 2 x e

x=-10:.1:10;

y=sin(3*x).*exp(-cos(x))+cos(2*x).*exp(-sin(x));

plot(x,y)

เราจะหาคําตอบในช่วง [0, 2]

x =0.95774795776341

MATLAB: Bisection Mtd.

function [x]=example91a(es)

% Calculate using Bisection Method between [0,2]

ea = inf;

xr = inf;

it=0;

xl=0;

xu=2;

while(ea > es)

it = it+1;

pxr=xr;

fxl=sin(3*xl).*exp(-cos(xl))+cos(2*xl).*exp(-sin(xl));

fxu=sin(3*xu).*exp(-cos(xu))+cos(2*xu).*exp(-sin(xu));

xr=(xl+xu)/2;

fxr=sin(3*xr).*exp(-cos(xr))+cos(2*xr).*exp(-sin(xr));

ea = abs((xr-pxr)/xr)*100;

x=[it xl fxl xu fxu xr fxr ea]

if(fxl*fxr > 0.0)

xl=xr;

elseif (fxl*fxr < 0.0)

xu=xr;

else

ea=0.0;

end

end

Bisection Results:>> example91a(0.01)

Iter xl fxl xu fxu xr fxr ea

x = 1.0000 0 1.0000 2.0000 -0.6869 1.0000 -0.0972 Inf

x = 2.0000 0 1.0000 1.0000 -0.0972 0.5000 0.7493 100.0000

x = 3.0000 0.5000 0.7493 1.0000 -0.0972 0.7500 0.4101 33.3333

x = 4.0000 0.7500 0.4101 1.0000 -0.0972 0.8750 0.1774 14.2857

x = 5.0000 0.8750 0.1774 1.0000 -0.0972 0.9375 0.0451 6.6667

x = 6.0000 0.9375 0.0451 1.0000 -0.0972 0.9688 -0.0249 3.2258

x = 7.0000 0.9375 0.0451 0.9688 -0.0249 0.9531 0.0104 1.6393

x = 8.0000 0.9531 0.0104 0.9688 -0.0249 0.9609 -0.0072 0.8130

x = 9.0000 0.9531 0.0104 0.9609 -0.0072 0.9570 0.0016 0.4082

x = 10.0000 0.9570 0.0016 0.9609 -0.0072 0.9590 -0.0028 0.2037

x = 11.0000 0.9570 0.0016 0.9590 -0.0028 0.9580 -0.0006 0.1019

x = 12.0000 0.9570 0.0016 0.9580 -0.0006 0.9575 0.0005 0.0510

x = 13.0000 0.9575 0.0005 0.9580 -0.0006 0.9578 -0.0000 0.0255

x = 14.0000 0.9575 0.0005 0.9578 -0.0000 0.9576 0.0002 0.0127

x = 15.0000 0.9576 0.0002 0.9578 -0.0000 0.9577 0.0001 0.0064

ans =

15.00000000000000 0.95776367187500 -0.00003535871565 0.95770263671875

0.00023929892750 0.95770263671875 0.00010197464576 0.00637308010962

x =0.95774795776341

True error = 0.004732%

Other Results: xt= 0.95774795776341

• Es = 0.01%

– It = 15;xr=0.95770263671875

– Ea=0.006373%, et = 0.004732%

• Es = Es = 0.001%

– It = 18;xr=0.95774078369141

– Ea=0.0007966%, et = 0.0007491%

• Es = 0.0001%

– It = 21;xr=0.95774745941162

– Ea=0.00009957%, et = 0.00005203%

• Es = 0.000001%

– It = 28;xr=0.95774795860052

– Ea=0.0000007779%, et = 8.740e-008%

MATLAB: False-Position

function [x]=example91b(es)

% Calculate using False-Position Method between [0,2]

ea = inf;

xr = inf;

it=0;

xl=0;

xu=2;

while(ea > es)

it = it+1;

pxr=xr;

fxl=sin(3*xl).*exp(-cos(xl))+cos(2*xl).*exp(-sin(xl));

fxu=sin(3*xu).*exp(-cos(xu))+cos(2*xu).*exp(-sin(xu));

% xr=(xl+xu)/2;

xr=xu-((fxu*(xl-xu))/(fxl-fxu));

fxr=sin(3*xr).*exp(-cos(xr))+cos(2*xr).*exp(-sin(xr));

ea = abs((xr-pxr)/xr)*100;

x=[it xl fxl xu fxu xr fxr ea]

if(fxl*fxr > 0.0)

xl=xr;

elseif (fxl*fxr < 0.0)

xu=xr;

else

ea=0.0;

end

end

FP Results:>> example91b(0.01)

Iter xl fxl xu fxu xr fxr ea

x = 1.0000 0 1.0000 2.0000 -0.6869 1.1856 -0.5611 Inf

x = 2.0000 0 1.0000 1.1856 -0.5611 0.7595 0.3940 56.1096

x = 3.0000 0.7595 0.3940 1.1856 -0.5611 0.9353 0.0500 18.7964

x = 4.0000 0.9353 0.0500 1.1856 -0.5611 0.9557 0.0045 2.1423

x = 5.0000 0.9557 0.0045 1.1856 -0.5611 0.9576 0.0004 0.1922

x = 6.0000 0.9576 0.0004 1.1856 -0.5611 0.9577 0.0000 0.0166

x = 7.0000 0.9577 0.0000 1.1856 -0.5611 0.9577 0.0000 0.0014

ans =

7.00000000000000 0.95773289766706 0.00003388653487 1.18559512875289

-0.56109590391892 0.95774665822935 0.00000292408716 0.00143676432376

x =0.95774795776341

True error = 0.0001357 %

Other Results: xt= 0.95774795776341

• สีนํ้าเงินได ้จาก Bisection Method

• สีเขียวได ้จาก False-Position Method

• Es = 0.01%

– It = 15;xr=0.95770263671875 Ea=0.006373%, et = 0.004732%

– It = 7;xr=0.95774665822935 Ea=0.001437%, et = 0.0001357%

• Es = Es = 0.001%

– It = 18;xr=0.95774078369141 Ea=0.0007966%, et = 0.0007491%

– It = 8;xr=0.95774784562942 Ea=0.0001240%, et = 0.00001171%

• Es = 0.0001%

– It = 21;xr=0.95774745941162 Ea=0.00009957%, et = 0.00005203%

– It = 9;xr=0.95774794808763 Ea=0.00001070%, et = 0.000001010%

• Es = 0.000001%

– It = 28;xr=0.95774795860052 Ea=0.0000007779%, et = 8.740e-008%

– It = 10;xr=0.95774795692851 Ea=0.0000009231%, et = 8.717e-008%

Newton-Ralphson Method

f ( x) = sin3

− cos

x e

x + cos2

− sin

x e

x

−cos

e

x

sin 3

−sin

x

e

x

x

x

x

cos

sin

cos2

f '( x) = sin 3 x d

+ e

d

+ cos2 x d

+ e

d

dx

dx

dx

dx

= sin3

−cos

x e

x ⋅ sin

−cos

x + e

x ⋅ 3cos3 x + cos2

−sin

x e

x (−cos x)

−sin

+ e

x ⋅ ( 2

− sin 2 x)

−cos

= e

x [sin 3 x sin x + 3cos3 x]

−sin

e

x [cos 2 x cos x + 2sin 2 x]

f ( x )

i

x

= x

i 1

+

i

f '( x )

i

sin 3

cos x

x e

i

+ cos2

sin x

x e

i

i

i

x

= x

i 1

+

i

−cos x

e

i [sin 3 x sin x + 3cos 3 x ]

sin x

e

i [cos 2 x cos x + 2sin 2 x ]

i

i

i

i

i

i

MATLAB Program:

function [x]=example91c(es,x0)

% Calculate solution using Newton-Ralphson, x0=initial;

it=0;

xi=x0;

ea=inf;

while (ea > es)

it = it+1;

fxi=sin(3*xi)*exp(-cos(xi))+cos(2*xi)*exp(-sin(xi));

dfxi=exp(-cos(xi))*(sin(3*xi)*sin(xi)+3*cos(3*xi))...

-exp(-sin(xi))*(cos(2*xi)*cos(xi)+2*sin(2*xi));

pxi=xi;

xi=pxi-fxi/dfxi;

ea=abs((xi-pxi)/xi)*100;

x=[it pxi fxi dfxi xi ea]

end

Result: ea=0.01, xo=?

• X0=0 โปรแกรมจะ Converge เข ้าสู่จุดอื่นด ้านซ ้าย

• X0=2 โปรแกรมจะ Converge เข ้าสู่จุดอื่นด ้านขวา

• ดูรูป

• ถ ้า x0 = 0.5 หรือ 1.5 โปรแกรมจะ Converge เข ้า

จุดที่ต ้องการอย่างรวดเร็วมาก

• เป็นไปได ้ที่เราเลือกจุดที่โปรแกรมไม่ Converge

• เราอาจจะใช ้ Bisection Method ก่อนเพื่อหาจุด

x0 ที่ดี จากนั้นต่อด ้วย Newton-Ralphson เพื่อให ้

ได ้คําตอบอย่างรวดเร็ว

เราจะหาคําตอบในช่วง [0, 2]

X0=0

x =0.95774795776341

2.1310

X0=2

Ex: หารากของ

− cos x

− sin x

f ( x) = sin 3 x e

+ cos 2 x e

x=-10:.1:10;

y=sin(3*x).*exp(-cos(x))+cos(2*x).*exp(-sin(x));

plot(x,y)

กรณีเลือก x0 = 0

-9.2837

Result: x0=0.5, es = 0.01

Iter xi fxi dfxi xi+1 ea

• x = 1.0000 0.5000 0.7493 -1.0485 1.2146 58.8351

• x = 2.0000 1.2146 -0.6362 -2.5824 0.9683 25.4413

• x = 3.0000 0.9683 -0.0238 -2.2755 0.9578 1.0939

• x = 4.0000 0.9578 -0.0001 -2.2502 0.9577 0.0061

• x = 4 0.95780651884294 -0.00013177280 -2.25024858111352

0.95774795962033 0.00611426231998

x =0.95774795776341

True error = 0.0000001939 %

Result: x0=0.5, es = 0.000001

Iter xi fxi dfxi xi+1 ea

• x = 1.0000 0.5000 0.7493 -1.0485 1.2146 58.8351

• x = 2.0000 1.2146 -0.6362 -2.5824 0.9683 25.4413

• x = 3.0000 0.9683 -0.0238 -2.2755 0.9578 1.0939

• x = 4.0000 0.9578 -0.0001 -2.2502 0.9577 0.0061

• x = 5.0000 0.9577 -0.0000 -2.2501 0.9577 0.0000

• x = 5 0.95774795962033 -0.00000000417826 -

2.25010587635298 0.95774795776341 0.00000019388349

x =0.95774795776341

True error = 0.000000000000000 %

Compare : xt =0.95774795776341

• Es = 0.01%

– It = 15;xr=0.95770263671875 Ea=0.006373%, et = 0.004732%

– It = 7;xr=0.95774665822935 Ea=0.001437%, et = 0.0001357%

– It = 4;xi= 0.95774795962033 Ea=0.006114%, et = 0.0000001939 %

• Es = Es = 0.001%

– It = 18;xr=0.95774078369141 Ea=0.0007966%, et = 0.0007491%

– It = 8;xr=0.95774784562942 Ea=0.0001240%, et = 0.00001171%

– It = 5;xi= 0.95774795776341 Ea=0.0000001939%, et < 1.0e-15 %

• Es = 0.0001%

– It = 21;xr=0.95774745941162 Ea=0.00009957%, et = 0.00005203%

– It = 9;xr=0.95774794808763 Ea=0.00001070%, et = 0.000001010%

– It = 5;xi= 0.95774795776341 Ea=0.0000001939%, et < 1.0e-15 %

• Es = 0.000001%

– It = 28;xr=0.95774795860052 Ea=0.0000007779%, et = 8.740e-008%

– It = 10;xr=0.95774795692851 Ea=0.0000009231%, et = 8.717e-008%

– It = 5;xi= 0.95774795776341 Ea=0.0000001939%, et < 1.0e-15 %

เพียง 5 iteration วิธีของ Newton-Ralphson

ให ้ Error น ้อยจน Double Precision วัดไม่ได ้

แต่ข ้อเสียคือจุด x0 จะต ้องเลือกให ้ดี

Chapter 10: System of

Linear Eq.

• จะ Limit อยู่ที่สมการ AX=B โดย A เป็น Square

Matrix

– N สมการ N Unknown

– จะมีคําตอบที่ Unique

– คําตอบจะมีได ้ต่อเมื่อ A ไม่เป็น Singular

• Determinant ไม่เท่ากับ 0

• A หา Inverse ได ้ และ X = A-1B

– ในกรณีที่ Determinant A ใกล ้ศูนย์ แต่ไม่ใช่ศูนย์

คําตอบจะ Sensitive กับ Error การคํานวณเมื่อมีการ

ปัดเศษจะต ้องระวัง

• กรณีนี้ เราเรียกว่ามันเป็น ‘Ill-Conditioned System’

System of Linear Equations

Krammer’s Rule

Solution ของ AX=C

• A-1AX=A-1C

• X=A-1C

• Inverse หาได ้ยาก แม ้จะใช ้ Computer

คํานวณ เพราะเป็น O(n4)

Solution by Elimination

Gauss Elimination

1. ใน Elimination Step จาก AX=C เราพยายามทําให ้

Matrix A อยู่ในรูป Upper Diagonal Matrix ด ้วย

ขบวนการ Elimination คือการบวกและลบแต่ละแถวเข ้า

ด ้วยกัน และค่า C จะถูกบวกลบตามไปด ้วย

2. เมื่อ A เป็น Upper Diagonal แล ้ว การแก ้สมการสามารถ

ทําได ้ง่าย โดยเราหา xn ก่อนในแถวสุดท ้ายของสมการ

จากนั้นนํา xn ที่หาได ้มาแทนค่า เพื่อหา xn-1 ในแถวรอง

สุดท ้าย เนื่องจากเป็นการแทนค่าเพื่อหา Unknown

ย ้อนหลัง เราจึงเรียก Back-Substitution

Gauss Elimination

• จาก

AX = C

Gauss Elimination

AX = C

a

a

a

x

c

1

,

1

,

1 2

,

1 n   1 

 1 

    

a

a

a

x

c

2 1

,

2,2

2, n   2 

 2

=

 

    

  

    

a

a

a

x

c

n 1

,

n,2

n, n   n

n

aˆ

aˆ

aˆ

x

cˆ

1

,

1

,

1 2

,

1 n   1 

 1 

    

a

aˆ

aˆ

x

cˆ

2 1

,

2,2

2, n   2 

 2

=

 

    

  

    

 0

0

aˆ

x

cˆ

n, n   n

n

Gauss Elimination

a

a

a

a

a x

c

1

,

1

,

1 2

,

1 3

,

1 4

,

1 5

 1  1

    

a

a

a

a

a

x

c

2 1

,

2,2

2,3

2,4

2,5   2 

 2

a

a

a

a

a   x  =  c

 3 1,

3,2

3,3

3,4

3,5   3  3

a

a

a

a

a

x

c

4 1

,

4,2

4,3

4,4

4,5   4 

 4

 

a

a

a

a

a

x

c

 5 1,

5,2

5,3

5,4

5,5   5 

 5

a

a

a

a

a x

c

1

,

1

,

1 2

,

1 3

,

1 4

,

1 5

 1  1 

/

/

/

/

    

 0

a

a

a

a

x

c'

2,2

2,3

2,4

2,5   2 

 2

/

/

/

/

0

a

a

a

a   x  =  c' 

3,2

3,3

3,4

3,5 

/

/

/

/

 3  3

 0

a

a

a

a

x

c'

4,2

4,3

4,4

4,5   4 

 4

/

/

/

/

 

a

a

a

a

x

0

c' 

5,2

5,3

5,4

5,5   5 

 5 

Gauss Elimination

a

a

a

a

a x

c

1

,

1

,

1 2

,

1 3

,

1 4

,

1 5

 1  1 

/

/

/

/

    

 0

a

a

a

a

x

c'

2,2

2,3

2,4

2,5   2 

 2

/

/

/

/

0

a

a

a

a   x  =  c' 

3,2

3,3

3,4

3,5 

/

/

/

/

 3  3

 0

a

a

a

a

x

c'

4,2

4,3

4,4

4,5   4 

 4

/

/

/

/

 

a

a

a

a

x

0

c' 

5,2

5,3

5,4

5,5   5 

 5 

a

a

a

a

a x

c

1

,

1

,

1 2

,

1 3

,

1 4

,

1 5

 1  1 

/

/

/

/

    

 0

a

a

a

a

x

c'

2,2

2,3

2,4

2,5   2 

 2 

//

//

//

0

0

a

a

a   x  =  c'' 

3,3

3,4

3,5 

//

//

//

 3  3

 0

0

a

a

a

x

c''

4,3

4,4

4,5   4 

 4

//

//

//  

a

a

a

x

0

0

c'' 

5,3

5,4

5,5   5 

 5 

Gauss Elimination

a

a

a

a

a x

c

1

,

1

,

1 2

,

1 3

,

1 4

,

1 5

 1   1 

/

/

/

/

   

 0

a

a

a

a

x

c'

2,2

2,3

2,4

2,5   2 

 2 

//

//

//

0

0

a

a

a   x  =  c'' 

3,3

3,4

3,5 

//

//

//

 3  3

 0

0

a

a

a

x

c''

4,3

4,4

4,5   4 

 4

//

//

//  

a

a

a

x

0

0

c'' 

5,3

5,4

5,5   5 

 5 

a

a

a

a

a x

c

1

,

1

,

1 2

,

1 3

,

1 4

,

1 5

 1   1 

/

/

/

/

   

 0

a

a

a

a

x

c'

2,2

2,3

2,4

2,5   2 

 2 

//

//

//

0

0

a

a

a   x  =  c'' 

3,3

3,4

3,5 

///

///

 3  3 

 0

0

0

a

a

x

c'''

4,4

4,5   4 

 4

///

///  

a

a

x

0

0

0

c''' 

5,4

5,5   5 

 5 

Gauss Elimination

a

a

a

a

a x

c

1

,

1

,

1 2

,

1 3

,

1 4

,

1 5

 1   1 

/

/

/

/

   

 0

a

a

a

a

x

c'

2,2

2,3

2,4

2,5   2 

 2 

//

//

//

0

0

a

a

a   x  =  c'' 

3,3

3,4

3,5 

///

///

 3  3 

 0

0

0

a

a

x

c'''

4,4

4,5   4 

 4

///

///  

a

a

x

0

0

0

c''' 

5,4

5,5   5 

 5 

a

a

a

a

a x

c

1

,

1

,

1 2

,

1 3

,

1 4

,

1 5

 1   1 

/

/

/

/

   

 0

a

a

a

a

x

c'

2,2

2,3

2,4

2,5   2 

 2 

//

//

//

0

0

a

a

a

  x  =  c'' 

3,3

3,4

3,5 

(3)

(3)

 3  3

(3) 

 0

0

0

a

a

x

c

4,4

4,5   4 

 4 

(4)  

a

x

 (4)

c

 0

0

0

0

5,5   5 

 5 

Back Substitution

a

a

a

a

a   x

c

1

,

1

,

1 2

,

1 3

,

1 4

,

1 5

1

1

   

/

/

/

/

0

a

a

a

a

x

c'

2,2

2,3

2,4

2,5 

2

2

  

//

//

//

 0

0

a

a

a   x  =  c'' 

3,3

3,4

3,5

3

3

   

(3)

(3)

(3)

 0

0

0

a

a

x

c

4,4

4,5

 4  4 

(4)

(4)

0

0

0

0

a

x   c

5,5   5 

 5 

a x + a x + a x + a x + a x = c 1

,

1

1

,

1 2

2

,

1 3 3

,

1 4

4

,

1 5 5

1

/

/

/

/

a x + a x + a x + a x = c'

2,2

2

2,3 3

2,4

4

2,5 5

2

//

//

//

a x + a x + a x = c''

3,3 3

3,4

4

3,5 5

3

(3)

(3)

(3)

a x + a x = c

4,4

4

4,5 5

4

(4)

(4)

a

x = c

5,5

5

5

Back

a x + a x + a x + a x + a x = c 1

,

1

1

,

1 2

2

,

1 3 3

,

1 4

4

,

1 5 5

1

/

/

/

/

+

+

+

=

Substitution

a x

a x

a x

a x

c'

2,2

2

2,3 3

2,4

4

2,5 5

2

//

//

//

a x + a x + a x = c''

3,3 3

3,4

4

3,5 5

3

(3)

(3)

(3)

a x + a x = c

4,4

4

4,5 5

4

(4)

(4)

a

x = c

5,5

5

5

(4)

c 5

x =

5

(4)

a 5,5

(3)

(3)

c

a x

4

4,5 5

x =

4

(3)

a 4,4

//

//

c'' − a x a x

3

3,4

4

3,5 5

x =

3

//

a 3,3

/

/

/

c' − a x a x a x

2

2,3 3

2,4

4

2,5 5

x =

2

/

a 2,2

c a x a x a x a x

1

,

1 2

2

,

1 3 3

,

1 4

4

,

1 5 5

x =

1

a 1,

1

Gauss Elimination

Gauss Elimination

Gauss Elimination

Gauss Elimination

Gauss Elimination

Gauss Elimination Alg

• Elimination by Forward Substitution

Gauss Elimination Alg

• Back-Substitution

Gauss Elimination Alg

• Back-Substitution

Example 8.1

Example 8.1

Example 8.1

• Back Substitution

ปัญหาของ Gauss Elimination

ปัญหาของ Gauss Elimination

ปัญหาของ Gauss Elimination

ปัญหาของ Gauss Elimination

Gauss-Jordan Method

Gauss-Jordan Method

Example 8.2

Example 8.2

 3

− 0 1

.

− 2

.

0

|

7 85

.

 .

0 1

7

− 3

.

0

|

3

.

19

 .

0 3

− 0 2

.

10

|

4

.

71

( R )

2 − ( R )

1 × .1/ 3 → ( R )

2

 3

− .

0 1

− .

0 2

|

85

.

7

 0

0033

.

7

− 0 2933

.

|

5617

.

19

 3

.

0

− 0 2

.

10

|

.

71 4

( R )

3 − ( R )

1 × 3

. / 3 → ( R )

3

3

− 0 1

.

− 2

.

0

|

85

.

7

0

7 0033

.

− 0 2933

.

|

5617

.

19

0 − 0 1900

.

10 0200

.

|

6150

.

70

( R )

1 = ( R )

1 / 3

1 − 0333

.

0

− .

0 0667 |

6167

.

2

0

0033

.

7

− .

0 2933 |

5617

.

19

0 − 0 1900

.

.

10 0200

|

6150

.

70

Example 8.2

1.0000 − 0.0333 − 0 0667

.

|

.

2 6167 

 0

7 0033

.

− 0 2933

.

|

− .

19 5617

 0

− .

0 1900

.

10 0200

|

.

70 6150 

( R )

3 = ( R )

3 − ( R )

2 × (−.

)

19 / 7 0033

.

 .

1 0000

− .

0 0333

− 0 0667

.

|

.

2 6167 

 0

7.0033

− 0.2933 | −

5617

.

19

 0

0

.

10 0120

|

.

70 0843 

( R )

2 = ( R )

2 / 7 0033

.

1.0000 − .

0 0333

− .

0 0667 |

.

2 6167 

 0

1

− 0.0419 | − .

2 7932

 0

0

.

10 0120

|

0843

.

70

( R )

1 = ( R )

1 − ( R )

2 × (−

)

0333

.

 .

1 0000

0

− 0 0681

.

|

2 5236

.

 0

1.0000

− .

0 0419 |

− 2.7932

 0

0

0120

.

10

|

.

70 0843 

Example 8.2

 0000

.

1

0

− 0 0681

.

|

5236

.

2

 0

0000

.

1

− 0419

.

0

|

− 7932

.

2

 0

0

0120

.

10

|

70 0843

.

( R )

3 = ( R )

3 /

012

.

10

 0000

.

1

0

− 0681

.

0

|

5236

.

2

 0

0000

.

1

− 0419

.

0

|

− 7932

.

2

 0

0

1 0000

.

|

0000

.

7

( R )

2 = ( R )

2 − ( R )

3 * (−

)

0419

.

 0000

.

1

0

− 0681

.

0

|

5236

.

2

 0

.

1 0000

0

|

− 5000

.

2

 0

0

0000

.

1

|

0000

.

7

( R )

1 = ( R )

1 − ( R )

3 * (−

)

0681

.

 0000

.

1

0

0

|

0000

.

3

 0

0000

.

1

0

|

− 5000

.

2

 0

0

0000

.

1

|

0000

.

7

Example 8.2

 .

1 0000

0

0

|

.

3 0000 

0

.

1 0000

0

|

− 2 5000

.

 0

0

.

1 0000 |

.

7 0000 

A'X = C'

1

0

0  x

 3 

1

   

0

1

0

x

= − .

2 5

2

   

0 0 1  x   7 

  3 

x

= 3

1

x

= 2

− 5

.

2

x

= 7

3

Example 8.2

การหา Matrix Inverse ด ้วย GJ

การหา Matrix Inverse ด ้วย GJ

การหา Matrix Inverse ด ้วย GJ

การหา Matrix Inverse ด ้วย GJ

การหา Matrix Inverse ด ้วย GJ

Example 8.3

Example 8.3

Example 8.3

1

3

− 1

.

− 2

. −

 3325

.

.0049

0068

.

1

.

7

− 3

.

=

− .0052 .1429 .0042

.3 − 2

.

10 

− 0101

.

.0027 .0999

จํานวนการคํานวณจะใช ้น ้อยกว่าวิธีทาง Analytic Method มาก

AdjA

A 1 =

A

Iterative Method and

Gauss-Seidel

Gauss-Seidel

Gauss-Seidel

Gauss-Seidel

Gauss-Seidel

Gauss-Seidel

Gauss-Seidel: Ex 8.4

Gauss-Seidel: Ex 8.4

Gauss-Seidel: Ex 8.4

0

0

x = x = x =

1

2

0 0

3

+

x +

x

+

+

1

.

7 85

1

.

0

0

2

.

0

0

85

.

7

(

1

.

0

)

0

.

0 (

2 )

0

.

7 85

2

3

x =

=

=

= .

2 61667

1

3

3

3

x +

x

+

1

.

19 3

.

0 1 1

3

.

0

0

19 3

.

0 (

1

.

.

2 61667)

.

0 (

3 0)

.

19 56167

1

3

x =

=

=

= −2 79452

.

2

7

7

7

x +

x

+

1

.

71 4

3

.

0

1

0 2

.

1

.

71 4 0 (

3

. 2 61667

.

)

(

2

.

0

)

79452

.

2

056095

.

70

1

2

x =

=

=

= 00561

.

7

3

10

10

10

2 61667

.

− 0

e

=

×100 =

%

100

a 1

,

61667

.

2

− 2 79452

.

− 0

e

=

×100 =

%

100

a 1

,

− 79452

.

2

7 00561

.

− 0

e

=

×100 =

%

100

a 1

,

00561

.

7

Gauss-Seidel: Ex 8.4

1

x = 2

x = −

x =

1

61667

.

, 1

.

2

,

79452

1

.

7 00561

2

3

+

x +

x

+

+

2

7 85

.

.

0 1 1

.

0 2 1

7 85

.

.

0 (

1

.

2

)

79452

0 (

2

.

)

00561

.

7

8 971696

.

8

2

3

x =

=

=

= 99056

.

2

1

3

3

3

x +

x

+

2

.

19 3 0 1

.

2

3

.

0

1

3

.

19

0. (

1

)

99056

.

2

(

3

.

0

7

)

00561

.

17 49737

.

1

3

x =

=

=

= − 49962

.

2

2

7

7

7

x +

x

+

2

71 4

.

0 3

.

2

.

0 2 2

71 4

.

(

3

.

0

)

99056

.

2

(

2

.

0

2

)

49962

.

70 00291

.

1

2

x =

=

=

= 00029

.

7

3

10

10

10

.

2 99056 − .

2 61667

e

=

×100 =

%

50

.

12

a,2

99056

.

2

− 49962

.

2

− (− .

2

)

79452

e

=

×100 =

%

80

.

11

a,2

− 49962

.

2

00029

.

7

− 7 00561

.

e

=

×100 =

%

076

.

0

a,2

7 00029

.

Gauss-Seidel: Ex 8.4

2

x = 2.

x = −

x =

1

,

99056

2

,

49962

.

2

2

00029

.

7

2

3

+

x +

x

+

+

3

7.85 0 1

.

2

0 2

.

2

85

.

7

(

1

.

0

)

49962

.

2

.

0 (

2

)

00029

.

7

0001

.

9

2

3

x =

=

=

= 3 00003

.

1

3

3

3

x +

x

+

3

3

.

19

1

.

0

3

3

.

0

2

3

.

19

)

00003

.

3

(

1

.

0

.

0 (

3 .

7

)

00029

499916

.

17

1

3

x =

=

=

= − .

2 499988

2

7

7

7

x +

x

+

3

.

71 4 0 3

.

3

.

0 2 3

4

.

71

.

3

(

3

.

0

)

00003

(

2

.

0

2

)

499988

.

.

69 99999

1

2

x =

=

=

= .

6 999999

3

10

10

10

00003

.

3

− 99056

.

2

e

=

×100 =

%

3157

.

0

a,3

00003

.

3

− 499988

.

2

− (−

)

49962

.

2

e

=

×100 =

%

01472

.

0

a,3

− 499988

.

2

999999

.

6

− 00029

.

7

e

=

×100 =

%

004157

.

0

a,3

999999

.

6

Gauss-Seidel: Ex 8.4

Gauss-Seidel: Ex 8.4

Jacobi Method

Convergence of Iterative

Method

Break

• After Break

– LU Decomposition

LU Decomposition

LU Decomposition

LU Decomposition

LU Decomposition

LU Decomposition

LU Decomposition

LU Decomposition

LU Decomposition

LU Decomposition

LU Decomposition

Crout Decomposition

Crout Decomposition

Crout Decomposition

Crout Decomposition

Crout Decomposition

Crout Decomposition

Crout Decomposition

Example 8.6

Crout Decomposition

Crout Decomposition

Example 8.7

Example 8.7

Summary Chapter 8

Homework 9, Chapter 10

DOWNLOAD

คํานวณแนะนําให้เขียนโปรแกรม หรือใช้

MATLAB หรือใช้ Spreadsheet(Excel)

หยุด 2 ส ัปดาห์

การบ้าน ส่งต้นช่ัวโมง พุธ 20 เมษายน

End of Chapter 10

• Next Week (Wk 13)

– No Class วันจักรี

• Following Week (Wk 14)

– No Class วันสงกรานต์

• Week 15; Chapter 11

– April 20

– Numerical Differentiation and

Integration

CPE 332

Computer Engineering

Mathematics II

Part III, Chapter 11

Numerical Differentiation and Integration

Today Topics

• Chapter 11 Numerical Differentiation and

Integration

– Derivative Approximation

• Forward, Backward, Centered Difference

• High Order Derivative

• High Accuracy Approximation

– Integral Approximation

• Polynomial

– Zero Order

– First Order (Trapezoidal)

– Second Order (Simpson 1/3)

– Third Order (Simpson 3/8)

• More Accurate Method

– Richardson Extrapolation

Slope of Line

m=slope = tan θ

= (y2-y1)/ (x2-x1)

(x2,y2)

(y2-y1)

θ

θ

(x1,y1)

(x2-x1)

Definition of Derivative

• Derivative of f(x) at any point x is the slope

of the tangent line at that point (x,f(x))

• Mathematically

y

f ( x + x

∆ ) − f ( x

Derivative of

( ) ≡ lim

=

)

f x

lim

x

∆ →0

x

x

∆ →0

x

• For function y=f(x), derivative of function is

written in many forms

dy

' ( ),

',

(read '

dee -

d

y ee - x' or

)

f

x

y

y

dx

Approximation of slope at point x

using secant line

• Slope at ‘ x’ ≅ [ f( x+∆ x) - f( x)] / ∆ x = y / x As ∆ x approaches 0,

we have ∆y/∆x

approach a true slope

of f at x.

( x+x, f( x+x))

f( x+x) - f( x)

( x, f( x))

=∆ y

x

Derivative(Forward)

f ( x + ∆) − f x

f '( x) =

( )

lim

∆→0

y

Derivative(Backward)

f ( x) − f ( x − ∆

f '( x) =

)

lim

∆→0

y

Derivative(Central)

f ( x + ∆) − f ( x − ∆

f '( x) =

)

lim

∆→0

2∆

y

x

x

Introduction

• ในทางปฎิบัติ เราได ้ Sample ของ Data เป็นจุด การหา

Derivative ก็คือการลบค่าของจุด Data ที่อยู่ติดกันและ

หารด ้วยระยะห่างระหว่างจุด นี่คือวิธีการของ Finite

Divided-Difference h h

f(xi+1)

f(xi)

f(xi-1)

xi-1 x

i xi+1

Finite Divided-Difference

Finite Divided-Difference

Finite Divided-Difference

f '''( x ) 3

h

f ( x ) − f ( x ) = 2 f '( x ) h + 2

i

+...

i 1

+

i 1

i

!

3

Second Derivative

Second Derivative

High Accuracy Finite

Divided-Difference

High Accuracy Finite

Divided-Difference

Summary

Summary

Summary

Numerical Integration

• Newton-Cotes Integration Formula

• Zero-Order Approximation

• First-Order Approximation

– Trapezoidal Rule

• Second-Order Approximation

– Simpson 1/3 rule

• Third-Order Approximation

– Simpson 3/8 rule

• Romberg Integration

– Richardson Extrapolation

– Romberg Integration Algorithm

Integral Approximation

Integral Approximation

x0 x1 x2 x3 x4

xn-1 xn

Zero-Order Approximation

x0 x1 x2 x3 x4

xn-1 xn

Zero-Order Approximation

• a0=f(a)

x0 x1 x2 x3 x4

xn-1 xn

b

b

I =

f ( x) dx a dx = a [ b a] = f ( a)[ b a]

∫ 0

0

a

a

Zero-Order Approximation

• a0=f((a+b)/2)

x0 x1 x2 x3 x4

xn-1 xn

b

b

a + b

I =

f ( x) dx a dx = a [ b a] = f (

)[ b a]

∫ 0

0

2

a

a

Zero-Order

App

roximation

First-Order Approximation

x0 x1 x2 x3 x4

xn-1 xn

First-Order Approximation

First-Order Approximation

First-Order Approximation

Trapezoidal Rule

Trapezoidal Rule

Second-Order Approximation

x0 x1 x2 x3 x4

xn-1 xn

Simpson’s 1/3 Rule

Second-Order Approximation

Second-Order Approximation

Second-Order Approximation

3n

Second-Order Approximation

Third-Order Approximation

x0 x1 x2 x3 x4

xn-1 xn

Simpson’s 3/8 Rule

Third Degree Approximation

Simson’s 3/8 Rule

Romberg Integration

Romberg Integration

Romberg Integration

Romberg Integration

Romberg Integration

Romberg Integration

Romberg Integration

Romberg Integration

Romberg Integration

Romberg Integration

Romberg Integration

Romberg Integration

Summary

Homework 10: Ch 11

• Download HW

• Next Week

– ส่งการบ ้าน HW 10

– Chapter 11 Solutions of ODE

– HW 11 (Option)

CPE 332

Computer Engineering

Mathematics II

Part III,

Chapter 12 ODE

Today Topics

Chapter 12 Ordinary Differential

Equation

HW 10 Due Today

HW 11 Due Wed 4 May, 12.00 Noon

ส่งห้องเลขาภาค

Course Ends This Week

No Chapter 13(Curve Fitting)

Differential Equation

• เป็นสมการทางคณิตศาสตร์ ที่ประกอบด ้วย

Derivative (หรือ Integral)

– Solution ของสมการคือ Function ที่ไม่มี Derivative

หรือ Integral ปรากฏอยู่

• รูปแบบทั่วไปคือ (First Order)

– dy/dx = f(x,y)

– กรณีพิเศษที่ dy/dx = f(x) เราสามารถแก ้สมการโดย

ใช ้การ Integrate

• dy = f(x)dx

• y = ∫ f(x)dx = F(x) + C : ค่า C สามารถหาได ้โดยการกําหนด

Initial condition

– ปกติการแก ้สมการทั่วไปไม่สามารถใช ้ Integrate ได ้

Example 1

• Solve for dy = 3 2 x; x = ,2 y = 4

dx

Example 1

• Solve for dy = 3 2 x; x = ,2 y = 4

dx

dy = 3 2

x dx

y = 3 2

x dx

3

y = x + C (

General So

lution)

Example 1

• Solve for dy = 3 2 x; x = ,2 y = 4

dx

dy = 3 2

x dx

y = 3 2

x dx

3

y = x + C (

General So

lution)

x = 2, y = 4

4 = 23 + C

C = 4

3

y = x − 4

Differential Equation

• ค่า Order สูงสุดของ Derivative จะ

กําหนด Order ของ Differential

Equation

– จํานวน Initial Condition ที่จะต ้องใช ้จะ

เท่ากับค่าของ Order ของสมการ

Example 2

d

F ( x, y) :

3

xy + 2

2

y + x = 4

Fin

d f ( x, y) =

F ( x, y)

dx

d

d

3

xy + 2

2

y + x =

4

dx

dx

d

dx

dy

dx

[

3

3

x

y + y

] + 2

+ 2 x

= 0

dx

dx

dx

dx

dy

dy

3

[

2

3

xy

+ y ]+ 2

+ 2 x = 0

dx

dx

dy

dy

3

2

xy

+ 2

= 2

3

x y

dx

dx

dy 3

[

2

xy + 2] = 2

3

x y

dx

dy

− 2

3

x

=

y

dx

3

[

2

xy + ]

2

Example 2

3

• แก ้สมการ dy

− 2 x

=

y

f ( x, y) =

dx

3

[

2

xy + 2]

− 2 x − 3

dy =

y dx

3

[ xy 2 + 2]

2

3

y = ∫ − x y dx

3

[ xy 2 + ]

2

Integrate ไม่ได ้

Differential Equation

• Tools ทางคณิตศาสตร์ที่สําคัญสําหรับใช ้ใน

การแก ้ปัญหา Differential Equation

– Laplace Transform (one side/two side)

– Fourier Transform ใช ้ได ้เช่นกัน

– Z-Transform ใช ้สําหรับแก ้ปัญหา Discrete

Version ของ Differential Equation

• คือ Differential Equation ที่ได ้จากการสุ่ม

ตัวอย่างของตัวแปร

• กรณีนี้เราเรียก Difference Equation

– ทั้งหมดนี้อยู่ในเนื้อหาวิชา CPE 308

Differential Equation

• บทนี้เราจะมาดู Numerical Method

สําหรับใช ้ในการแก ้สมการ Differential

Equation

– เราจะจํากัดอยู่ที่ First Order และ Initial

Condition ที่จุดตั้งต ้น

– เป็นสมการของ Ordinary Differential

Equation

• One Independent Variable

– สามารถดัดแปลงสําหรับ Higher Order ได ้

– กรณีที่เป็น Boundary Condition จะต ้องใช ้วิธี

อื่น

• ศึกษาเพิ่มเติมได ้จาก Reference

Chapter 10

• Ordinary Differential Equation

ODE

ODE

ODE

ODE: One Step Method

y

y 4

3

y 1

y

y

2

0

x

h

0

x

x

x

x

x

0

1

2

3

4

ODE: One Step Method

Euler’s Method

Euler’s Method

Euler’s Method

Euler’s Method

Euler’s Method

Euler’s Method

Heun and Polygon Method

Heun Method

Heun Method

Heun Method

Heun and Polygon Method

Heun Method

Heun and Polygon Method

Heun Method

Heun Method

Improved Polygon Method

Improved Polygon Method

Runge-Kutta Method

Runge-Kutta Method

Second Order Runge-Kutta

Method

Second Order Runge-Kutta

Method

Second Order Runge-Kutta

Method

Second Order Runge-Kutta

Method

Second Order Runge-Kutta

Method

Second Order Runge-Kutta

Method

Second Order Runge-Kutta

Method

Second Order Runge-Kutta

Method

Second Order Runge-Kutta

Method

Third Order Runge-Kutta

Method

Third Order Runge-Kutta

Method

Forth Order Runge-Kutta

Method

Forth Order Runge-Kutta

Method

Higher Order Runge-Kutta

Method

Comparison

Chapter 12 Homework (HW 11)

Download HW 11(ODE)และทําใน

Sheet

Option ถ้าใครส่งครบ 10 ครั้งและได้เต็ม

ไม่ต้องส่ง

จะนับ 10 HW ที่คะแนนมากที่สุด

ส่งพุธที่ 4 พ.ค. ที่ห้องภาค 5-310 ก่อน

เที่ยง

ไม่มีบทที่ 13 Curve Fitting

Course Ends

Final Exam Preparation

สูตรจะให้มา

ข้อสอบมี 7 ข้อ เลือกทํา 5 ข้อ 10x5 =

50 คะแนน เทียบเป็นคะแนนเก็บ 50%

2 ข้อ เป็นเรื่องก่อน Midterm (Part 1

หนึ่งข้อ และ Part 2 หนึ่งข้อ)

5 ข้อเป็นเรื่องใหม่ ด ังน้ี

5 ข้อ เป็นเรื่องใหม่ หล ัง MT ด ังนี้

– 1. Taylor Series และการประมาณค่าของ

Function รวมถึง Error

– 2. Root of Function (Bisection หรือ

Newton)

– 3. Linear Equation 1 ข ้อ

• Gauss Elimination, Gauss Jordan, Gauss Seidel

และ LU Decompositon

– 4. Numerical Integration (Trapezoidal หรือ

Simpson)

– 5. ODE โดยใช ้ 4th Order RK

Formulas

ઙ્ગહ્ણહ્વ� MT+ઙ્ઘ�

�����

CPE 332

Computer Engineering

Mathematics II

Mathematics II

Chapter 12

Curve Fitting

(Last)

Curve Fitting

• บ่อยๆครั้ง เราได ้ข ้อมูลจากการวัดในสนาม ซึ่งมี

ลักษณะเป็นจุด และเราจะต ้องการหา Function

ทางคณิตศาสตร์เพื่อเชื่อมต่อจุดเหล่านี้

• Function ง่ายๆที่ใช ้กันคือ Polynomial ที่

Degree ต่างๆ

• Function ที่สลับซับซ ้อนขึ้นได ้แก่ Cosine

Function

• กรรมวิธีเหล่านี้ เพื่อที่จะหาค่าของ Function ที่อยู่

ระหว่างจุดที่เราวัด ซึ่งเราเรียก Interpolation

Interpolation

y

x

Linear Interpolation

y

x

Second Degree Piecewise

Interpolation

y

x

Third Degree Piecewise

Interpolation

y

x

Higher Order Degree

• ไม่นิยมใช ้ เพราะจะเกิดการ Oscillate ได ้

• Function ที่รอยต่อของ Polynomial จะไม่

ต่อเนื่อง

• N Degree จะใช ้ N+1 จุดสําหรับหาค่าของ

Function

Spline Interpolation

• เป็นการใช ้ Polynomial Degree ตํ่าๆ ที่นิยมคือ Third

Degree หรือ Cube

– เรียก Cubic Spline

• Cubic Spline จะหา Polynomial จากทีละ 2 จุด(ไม่ใช่ 4)

แต่บังคับให ้รอยต่อระหว่างแต่ละ Polynomial มีความ

ต่อเนื่อง โดยกําหนดให ้ Polynomial ที่ต่อกันมี First

Derivative และ Second Derivative เท่ากัน

• 2 จุดได ้ 2 สมการ บวกกับอีกสองสมการของ First และ

Second Derivative เป็น 4 สมการ 4 Unknown

• ที่ปลายและหัวจะสมมุติ First Derivative จากสมการ

เส ้นตรง และสมมุติ Second Derivative เท่ากับศูนย์ วิธีนี้

เราเรียก Natural Spline

• เป็นวิธีที่นิยมมากกว่า

Spline Interpolation

y

x

การคํานวณ Cubic Spline

Interpolation (1)

• พิจารณาจากจุดของตัวอย่าง ( x , y ); i

i

i

= ,...,

0

n −1

• เราต ้องการลาก Curve ผ่านจุดเหล่านี้ โดยใช ้

Third Degree Polynomial

2

3

y = f ( x) = a + a x + a x + a x 3

0

1

2

3

– มี 4 Unknown ที่ต ้องหา ดังนั้นต ้องสร ้าง 4 สมการ และ

ใช ้ 4 จุดของ Data ในการแก ้

– ผลคือ 3rd Degree Polynomial จะลากผ่านทีละ 4 จุด

ดังนั้นจะมีรอยต่อทุกๆ 4 จุด ที่ไม่ราบเรียบ เนื่องจากแต่

ละ 4 จุดจะใช ้ Polynomial คนละตัวกัน

– Polynomial ที่ต่อกันจะใช ้จุดร่วมกันที่รอยต่อ

การคํานวณ Cubic Spline

Interpolation (2)

• Example: กําหนด Data 7 จุดดังนี้

– Yi = [2 0 3 4 5 3 1]; Xi = [0 1 2 3 4 5 6]

• จงคํานวณ Third Degree Polynomial ที่

จะ Fit Data เหล่านี้

การคํานวณ Cubic Spline

Interpolation (3)

y0

y

1

เราต้องใช้ Cubic Polynomial 2 ต ัวแรก Fit 4 จุดแรก

ต ัวที่ 2 Fit 4 จุดถ ัดไป โดยจุดที่ 4 เป็นจุดร่วม

Cubic Spline Int.(4)

• Polynomial ตัวแรก คํานวณได ้ดังนี้

2 = a 0

0 = a + a + a + a

0

1

2

3

3 = a + 2 a + 4 a + 8 a

0

1

2

3

4 = a + 3 a + 9 a + 27 a

0

1

2

3

• แก ้สมการได ้

2

3

y = 2 − 6.8333 x + 6 x −1.1667 x

0

Cubic Spline Int.(5)

• Polynomial ตัวที่สอง คํานวณได ้ดังนี้

4 = a + 3 a + 9 a + 27 a

0

1

2

3

5 = a + 4 a +16 a + 64 a

0

1

2

3

3 = a + 5 a + 25 a +125 a

0

1

2

3

1 = a + 6 a + 36 a + 216 a

0

1

2

3

• แก ้สมการได ้

2

3

y = −47 + 35 x − 7.5 x + 0.5 x

1

Cubic Spline Int.(6)

Cubic Spline Int.(7)

Cubic Spline Int.(8)

Cubic Spline Int.(9)

• สิ่งที่เราต ้องการคือทําให ้รอยต่อราบเรียบ

– โดยการกําหนดให ้ค่า First Derivative ที่ปลายของ

Polynomial ตัวแรก เท่ากับค่า First Derivative ที่จุด

ต ้นของ Polynomial ตัวที่ 2

– อีกนัยหนึ่ง คือบังคับให ้ค่า First Derivative ที่จุดร่วม

ของสอง Polynomial ที่มาต่อกัน มีค่าเท่ากัน

– ในกรณีเราสามารถสร ้าง 4 สมการ เพื่อหา 4

Coefficient ได ้จากเพียงค่าของ Function ที่สองจุด

และค่าของ First Derivative ที่สองจุดนั้น (4 ค่า รวม)

• กล่าวคือ แต่ละช่วงของ Data (ระหว่างสองจุด) เราจะคํานวณ

Third Degree Polynomial หนึ่งตัว

Cubic Spline Int.(10)

• พิจารณาจาก Data n จุดของ xii=1,..,n คือ yi;

i=1,..,n เราต ้องการหา Cubic Polynomial,

Si(x) ทั้งหมด n-1 function เพื่อทําการ

Interpolate ค่าระหว่างจุด ในแต่ละช่วงของ

สองจุดที่ติดกัน

s ( x); x x < x

1

1

2

 s ( x); x x < x

S( x)

2

2

3

= 

s ( x);

n

xn− ≤ x < x

1

1

n

s ( x)

i

= a ( x

i

x )3

i

+ b ( x

i

x )2

i

+ c ( x

i

x ) i + d ; i

i

= ,

1 ,...,

2

n −1

Cubic Spline Int.(11)

• สมการ Polynomial ที่ใช ้เป็นแค่การเปลี่ยน

Variable โดยเลื่อน (Delay) สมการจากจุด

xi ให ้มาอยู่ที่ตําแหน่ง ศูนย์

• ค่า First และ Second Derivative ของ

แต่ละ Polynomial สามารถแสดงได ้ดังนี้

s ( x) = a ( x x )3 + b ( x x )2 + c ( x x ) + d i

i

i

i

i

i

i

i

s '( x) = a

3 ( x x )2 + b

2 ( x x ) + c

i

i

i

i

i

i

s ''( x) = 6 a ( x x ) + b

2

i

i

i

i

Cubic Spline Int.(12)

• ถ ้าจุด xi เป็นจุดร่วมระหว่างสอง

Polynomial กล่าวคือ

s ( x ) = s ( x )

i

i

i 1

i

• และที่จุด xi เราจะได ้ค่าของ Function คือ

y = s ( x ) = a ( x x )3 + b ( x x )2 + c ( x x ) + d i

i

i

i

i

i

i

i

i

i

i

i

i

y = d

i

i

• นอกจากนี้ เราได ้(สมมุติว่าแต่ละจุดสุ่ม

ตัวอย่างด ้วยระยะห่างเท่ากัน)

y = s ( x )

= a ( x

x )3

+ b ( x

x )2

+ c ( x

x )

+ d

i

i 1

i

i 1

i

i 1

i 1

i

i 1

i 1

i

i 1

i 1

y = d = a h 3

+ b h 2

+ c h

+ d ; h

= x x ; i

=

,...,

3

,

2

n

i

i

i 1

i 1

i 1

i 1

i

i 1

Cubic Spline Int.(13)

• เราสร ้างสมการจากการกําหนด First

Derivative ของจุดต่อให ้เท่ากัน

s '( x ) = s '( x )

ดังนั้น เราได

i

i

i 1

i

s '

i

( x) = 3 a ( x

i

x )2

i

+ 2 b ( x

i

x ) i + ci

s ' ( x )

i

i

= ci = s '( x )

i

i

= 3

2

aih + 2 bih + c ;

i

i =

,...,

3

,

2

n −1

1

1

1

1

• นอกจากนี้แล ้ว ค่า Second Derivative

จะต ้องต่อเนื่องกันตลอดช่วง

s ''( x )

i

i+

= s ''( x );

i+

i+

i = ,

1 ,...,

2

n −1

1

1

1

Cubic Spline Int.(14)

• แต่ค่า s ''( x) = 6 a ( x x ) + b 2

i

i

i

i

• ดังนั้น

s ''( x )

b

i

i

= 2 i

• และ

s ''( x )

2 b

i 1

+

i

=

1

+

i 1

+

s ''( x )

+

= 6 a ( x + − x ) + b

2

= 6 a h + b

2

i

i 1

i

i 1

i

i

i

i

• จาก

s ''( x )

i

i +

= s ''( x );

i +

i +

i = ,

1 ,...,

2

n −1

1

1

1

• เราได ้

b

2 + = 6 a h + b

2

i 1

i

i

Cubic Spline Int.(15)

• เพื่อให ้สมการดูง่าย เรากําหนด

M = s ''( x )

i

i

i

• ดังนั้น

b = M / 2

i

i

• และที่ได ้ก่อนหน้านี้ d = y

i

i

• ส่วน ai และ ci สามารถเขียนใหม่ได ้เป็น

M +1 − M

b

2 +1 = 6 a h + b

2  a

i

i

=

i

i

i

i

6 h

d = a h 3

2

−1

+ b h

−1

+ c h

−1 + d

i

i

i

i

i −1

3

2

y +1 − y

M +1 + 2 M

d +1 = a h + b h + c h + d c i

i

i

i

=

− 

h

i

i

i

i

i

i

h

6

Cubic Spline Int.(16)

• เราสรุปสูตรในการคํานวณดังนี้

M +1 − M

i

i

a =

i

6 h

Mi

b =

i

2

y +1 − y

M +1 + 2 M

i

i

i

i

c =

− 

h

i

h

6

d = y

i

i

• เมื่อเขียนในรูปของ Matrix เราได ้

Cubic Spline Int.(17)

• จาก

s '( x

i

i+ ) = 3

2

a h

i

+ 2 b h

i

+ ci = s '( x )

i+

i+

= c

1

1

1

i 1

+

Mi+ − Mi

Mi

yi+ − yi Mi+ + Mi

yi+ − yi

+

M i+ + Mi

1

2

2

2

3

+

h

h

h

h

+ 2

1

1

2

1

2

1

+

− 

 =

− 

6 h

 2 

h

6

h

6

h ( Mi

yi yi+ +

+

y

M

i

i+ + M i+

=

+

1

2 )

2

4

1

2

6

h

M i +

6  yi − 2 yi+ + yi

4

1

+2

M i+ + Mi+ = 

;

i = ,

1

,..,

3

,

2

n − 2

1

2

h

h

• สามารถเขียนเป็นสมการ Matrix ได ้ดังนี้

Cubic Spline Int.(18)

M 1 

1 4 1 0  0 0 0 0 M

2

2

y 1 − y 2 + y 3

0

1

4

1  0 0 0 0  M

2

3

y 2 − y 3 + y 4 

0 0 1 4  0 0 0 0 M

2

4 

y 3 − y 4 + y

6 

5

        

=

2 

0 0 0 0  4 1 0 0

h

M

y

2 y

y

n−3

n −4 −

n −3 +

n −2 

0 0 0 0  1 4 1 0 M

y

2 y

y

n −2

n−3 −

n −2 +

n −1 

0 0 0 0  0 1 4 1 M

y

2 y

y

n −1 

n −2 −

n −1 +

n

Mn

• ระบบประกอบไปด ้วย n-2 สมการและ n unknown

ดังนั้นเราต ้องการอีก สอง Condition เพื่อที่

สามารถจะแก ้สมการได ้

– สอง Condition สามารถเลือกได ้หลายแบบ ทําให ้เกิด

Variation ของ Cubic Spline Interpolation หลาย

วิธี

Cubic Spline Int.(19)

• ในที่นี้จะกล่าวถึง Spline สามแบบ ตาม

สอง Condition ที่กําหนดเพิ่มเติม

– Natural Spline

• โดยการกําหนด M1=Mn=0

– Parabolic Runout Spline

• โดยการกําหนด M1=M2 และ Mn=Mn-1

– Cubic Runout Spline

• โดยการใช ้ M1=2M2-M3 และ Mn=2Mn-1-Mn-2

• ในที่นี้จะกล่าวรายละเอียดเฉพาะ Natural

Spline

Cubic Spline Int.(20)

• Natural Spline

– เมื่อให ้ M1=Mn=0 Matrix จะลดรูปเหลือ n-2

สมการ และ n-2 Unknown ดังนี้

4 1 0  0 0 0 M

2

2 

y 1 − y 2 + y 3 



1

4

1  0 0 0

M

y

2

3

2 −

y 3 + y



4

0 1 4  0 0 0 M

y 2

4

3 −

y 4 + y

5



6 

         =

2 

h

0 0 0  4 1 0 M

y

2 y

y

n−3

n−4 −

n−3 +



n−2 

0 0 0  1 4 1 M

y

2 y

y

n−2 

n−3 −

n−2 +

n−1 



0 0 0  0 1 4 M

y

2 y

y

n−1 

n−2 −

n−1 +

n

– เมื่อเราแก ้สมการได ้ค่า Mi เราสามารถคํานวณ

ค่า Coefficient ของแต่ละ Polynomial ได ้

Cubic Spline Int.(21)

• สังเกตุว่า Coefficient Matrix ของระบบมีค่าใน

ส่วน Diagonal ที่เป็นค่าคงที่ ซึ่ง Matrix นี้เป็น

Matrix พิเศษ เรียก Toeplitz Matrix

– Toeplitz Matrix สามารถแก ้ปัญหาได ้ โดยใช ้ O(n2)

Computation Time แทนที่จะเป็น O(n3)

• Algorithm ที่ใช ้คือ Levinson Algorithm

– Toeplitz Matrix สามารถทํา LU Decomposition

โดยใช ้ O(n2) เช่นเดียวกัน

• รายละเอียดของคุณสมบัติพิเศษของ Matrix นี้

และ Levinson Algorithm อยู่นอกขอบเขตของ

วิชานี้ ผู ้สนใจสามารถศึกษาเพิ่มเติมได ้จากตํารา

ด ้าน Linear Algebra และ Numerical Method

Example: Natural Spline

• จากตัวอย่างเดิม เราจะใช ้ Spline มาทําการ Fit

Data ในกรณีนี้ n=7, h = 1, M0=M7=0

– Yi = [2 0 3 4 5 3 1]; Xi = [0 1 2 3 4 5 6]

• เราแก ้สมการ Matrix ดังนี้

4 1 0 0 0 M 2

2 − 2⋅0 + 

3

 5 



1

4

1

0

0

M

0

3

− 2⋅3+ 4



 2

0 1 4 1 0 M

4

= 63− 2⋅4 + 5 = 6 0 



0 0 1 4 1 M 5 

4 − 2⋅5 + 

3

− 3

0 0 0 1 4 M

5 2 3 1

 0 

 6 

 − ⋅ + 

• และได ้คําตอบคือ

T

M = [ ,

0

.

8

,

9923 − .

5

,

9692 2

,

8846

.

5

1

,

5692

.

,

3923

.

0]

T

M = [ ,

0

.

8

,

9923 5

,

9692

.

.

2

,

8846 − .

5

.

1

,

5692

,

3923

]

0

• ดังนั้น Coefficient ของ Cubic ทั้ง 6 จะ

เป็น

M i M

1

ai =

+

i = 1

[ .4987,−2.

1

,

4936 .

,

4756 −1.

1

,

4090 .

,

1603 −.

]

2321

6 h

b

M

i =

/ 2

i

= [ ,

0 .

4

,

4962 − .

2

.

1

,

9846

,

4423 − .

2

0

,

7846 .

]

6962

yi+ − yi Mi+ + 2 Mi

1

1

ci =

− 

h

h

6

= [−3.4987,0.

,

9974 2.

,

5090 0.9667,−0.

,

3756 −2.

]

4641

d

y

i =

i = [

,

3

,

0

,

2

]

3

,

5

,

4

Polynomial Coefficient Matrix

>> p=[a b c d]

p =

1.4987 0 -3.4987 2.0000

-2.4936 4.4962 0.9974 0

1.4756 -2.9846 2.5090 3.0000

-1.4090 1.4423 0.9667 4.0000

1.1603 -2.7846 -0.3756 5.0000

-0.2321 0.6962 -2.4641 3.0000

ต่อไปเราจะลอง Plot แต่ละ Polynomial ในแต่ละช่วง

ท ั้งหมด 6 Polynomial และ 6 ช่วง

Natural Spline Interpolation

Natural Spline Interpolation = R;

First Der = G; Second Der = B

Natural Spline Interpolation = R;

First D = G; Second D = B vs 6th D Poly = K

ปกติแนวโน ้มของ Polynomial ที่มี Degree สูงจะเกิด Overshoot

แต่ Spline จะให ้ Curve ที่ราบเรียบกว่า และใกล ้เคียงความจริง

y = 0.0347 6

x − 0 6375

.

5

x + .

4 4931 4

x

.

15 3125 3

x +

.

25 4722 2

x

.

16 05 x + 2

Regression

• บางครั้ง Data ที่เราเก็บได ้ กระจาย และการ

ทํา Interpolation ด ้วย Polynomial จะไม่

เหมาะสม

• การกระจายของ Data อาจจะเกิดจากความ

ไม่แน่นอนในการวัด หรือเป็นลักษณะ

Random

• อาจจะเกิดจาก Noise ในการวัด

– การ Fit ด ้วย Polynomial หรือแม ้แต่ Spline

จะให ้คําตอบที่ไม่ตรงกับความจริง

Regression: Data with Noise

Y

X

Polynomial Fit Noise Data

Y

X

Ordinary Simple Linear

Least-Squares Regression

• ในการนี้ การประมาณค่าโดยใช ้ Function ง่ายๆ เช่น

เส ้นตรง หรือ Exponential จะเหมาะสมกว่า

– โดยเฉพาะถ ้าเรารู ้พฤติกรรมของระบบอยู่บ ้าง และคาดหวังว่า Data

ที่ได ้ควรจะมีลักษณะอย่างไร

• การหาสมการที่จะ Fit ที่สุด โดยมีค่า Error ตํ่าสุด เรา

เรียกเป็นการทํา Regression

• กรณีของสมการเส ้นตรง (สําหรับตัวแปรเดียว)เราเรียก

Linear Regression

– ที่ถูกควรเรียก Simple Linear Regression ใช ้สําหรับตัวแปรเดียว

จะเป็นสมการเส ้นตรงในสองมิติ

• เนื่องจากมี Multiple Linear Regression ด ้วยสําหรับกรณีหลายตัว

แปร ซึ่งผลจะเป็น Linear Equation ในหลายมิติ

– สมการเส ้นตรงคือ y=ax+b (บางตําราใช ้ y = a + bx ซึ่งในกรณีนี้

ต ้องสลับค่า a และ b ในการคํานวณที่จะกล่าวต่อไป)

– Parameter 2 ตัวที่ต ้องหาคือ a = slope และ b = y-intercept

Ordinary Simple Linear

Least Square Regression

• สมการที่ Fit ที่สุดคือให ้ Error รวมตํ่าสุด

เรามักจะใช ้วิธีของ Least-Square เรียก

Ordinary Least Square Regression

(OLS) (หรือ Linear Least Square)

– ผลรวมของ Error ยกกําลังสองของแต่ละจุด

ของ Data เทียบกับเส ้นตรงที่ใช ้ มีค่าตํ่าสุด

• Linear Regression ใช ้กันมากในทางสถิติ

เพื่อหาความสัมพันธ์แบบเส ้นตรงของ

Random Variable สองตัว

• Linear Regression ประเภทอื่นๆก็มีใช ้กัน

อยู่ แต่ OLS จะนิยมที่สุด

Linear Least-Squares

Regression

Y

Error

X

Ordinary Linear Least

Square Regression (OLS)

• ค่า a และ b หาได ้จากสมการ Summation ของ

Error ยกกําลังสอง ของทุกจุด

– ค่าที่ตํ่าสุดหาได ้จากการหาค่า Derivative ของสมการ

และตั้งให ้เท่ากับ 0 (จุด Minima)

• ทั้งสอง Unknown คือ a และ b หาได ้จากสมการ

Derivative สองสมการที่ตั้งให ้เท่ากับ 0 ดังนี้

– df(a,b)/da = 0

– df(a,b)/db = 0

– f(a,b) คือ Sum of Square Error Function

• อีกค่าที่ใช ้ในทางสถิติคือค่า r = correlation

coefficient เป็นตัวบอกว่าเส ้นตรงที่ใช ้มัน Fit

ได ้ดีเพียงไร, r จะมีค่าระหว่าง [-1,1]

Linear Least-Squares

Regression

e

y

ax

b i

n

i =

i

i

; = ,

1 ,...,

2

Regression line

n

n

f ( a, b) =  2

e

y

ax

b

(x

i

= ( i i − 2)

i,yi)

i=1

i=1

yi

n

ei

df

axi+b

= −2[( y ax b x

i

i

) ]

i

= 0

da

i=1

n

df = −2( y ax b

i

i

) = 0

x

db

i

i=1

nx y

x

y

i

i −  i

a =

i

สมมุติเราได้ Data Pair (xi,yi);

n 2

x

x

i − (

2

)

i

i = 1,2,…,n ท ั้งหมด n จุด

b = y ax

Variable X และ Y อาจจะเป็น

nx y

x

y

i

i −  i

ได้ท ั้ง Continuous และ

r =

i

Discrete

n 2

x

x

n

y

y

i − (

2

)

i

 2 i −( 2)

i

n

n

− 2[( y ax b x

− 2[( y ax b)] = 0

i

i

) ]

i

= 0

i

i

i 1

=

i 1

=

n

n

n

n

n

2

y x a x b x

 −  − =

i

i

i −  i =

0

(1)

y

a

x

nb

0

i

i

i 1

=

i 1

=

i 1

=

i 1

=

i 1

=

r

eplace b in

(1)

1 n

1 n

b

n

n

n

n

n

=

y a x = Y aX

1

1

2

n

i

y x a x

y

a

x

x

i

i

i −  i −  

i

n

i

1

i

i = 0

=

i 1

=

n

n

i 1

=

i 1

=

i 1=

i 1

=

i 1=

2

n

1 n

n

n

1

n

2

y x

y

x

a

x

x

i

i

ii =

n

 i −  

i

n

i 1

=

i 1

=

i 1

=

i 1

=

i 1

=

 

2

n

n

n

n

n

2

ny x

y

x

a n

x

x

i

i −  i i =

  i −  

i

i 1

=

i 1

=

i 1

=

i 1

=

i 1

=

 

n

n

n

n

x y

x

y

i

i −  i i

i 1

=

i 1

=

i 1

a =

=

2

n

n

2

nx

x

i −  

i

i 1

=

i 1= 

R, r : Correlation Coefficient

ส่วน Correlation Coefficient คือ

Pearson Product-Moment Correlation Coefficient

นิยามว่าเท่าก ับอ ัตราส่วนของค่า Covariance ต่อผลคูณของค่า

Standard Deviation ของท ั้งสอง Variable (ของต ัวอย่าง)

และ สามารถแสดงได้ว่ามีค่า

r = CXY

S S

X

Y

nx y

x

y

i

i −  i

r =

i

n 2

x

x

n

y

y

i − (

2

)

i

 2 i − ( 2)

i

r = 1 หรือ -1 แสดงว่ามีความสัมพันธ์แบบเส ้นตรง

อย่างสมบูรณ์ (ในเชิงบวก หรือในเชิงลบ)

ถ ้า r เป็นศูนย์แสดงว่าไม่มีความสัมพันธ์กันเลย

Nonlinear Regression

• ในกรณีที่เส ้นตรงไม่เหมาะสมที่จะนํามาใช ้

เราสามารถจะใช ้ Function อื่นมาทํา

Regression เช่น

– Exponential Model

– Power Equation

– Saturation-Growth-Rate Equation

– Polynomial Regression

• รายละเอียดอยู่นอกเหนือวิชานี้

Example

• 1. จากการเก็บข ้อมูลตัวอย่าง วิชา CPE 332 ระหว่าง

คะแนน ที่นักศึกษาได ้ กับจํานวนรวมเป็นชั่วโมงที่

นักศึกษาขาดเรียน หรือมาสาย ได ้ข ้อมูลดังตาราง

ข ้างล่าง

Data

Hour

Grade

Data

Hour

Grade

1 5 72 6 10 40

2 8 51 7 7 68

3 2 86 8 6 63

4 6 75 9 12 46

5 4 88 10 8 65

– จงใช ้ OLS Regression แสดงสัมพันธ์ระหว่างตัวแปรทั้ง

สอง ทําการ Plot Scatter Diagram และเส ้นตรงที่ได ้

จากนั้นให ้หาค่า Correlation Coefficient

Example

Example

10

10

1

n =

;

10  x

X

x

i =

;

68

=

i = .68

i 1

10

=

i 1

=

10

1 10

y

Y

y

i =

;

654

=

i = .

65 4

i 1

10

=

i 1

=

10

10

10

x y

x

y

i

i =

;

4068

2

i = ;

538

2

i =

;

45084

i 1

=

i 1

=

i 1

=

n

n

n

n

x y

x

y

i

i −  i i

10

i=

i=

i=

× 4068 − 68×654

1

1

1

a =

=

= − 0159

.

5

2

n

n

10× 538 − 682

2

nx

x

i −  

i

i 1

=

i 1= 

b = Y a X =

.

65 4 − (− .

5

)

0159 × .

6 8 =

.

99 5079

nx y

x

y

i

i −  i

r =

i

2

nx

x

n

y

y

i − (

)2

2

i

i −( )2

i

10 × 4068 − 68× 654

=

= − 9069

.

0

10 × 538 − 682 10× 45084 − 6542

y=-5.0159x + 99.5079

r = -0.9069

End of Chapter 12

• Download

• Homework 12 ส่งก่อนวันพุธหน ้าที่ 3 ธ.ค.

ก่อน 12.00 น. (ส่งที่สาขา ใส่กล่อง ห ้อง 5-

310)

Course Ends

• Prepare For Exam

– ข ้อสอบมี 9 ข ้อ ทําทุกข ้อ (นําเครื่องคิดเลขมาด ้วย)

– เรื่องก่อน Midterm จะออก 3 ข ้อ จาก 6 บทแรก + 6

ข ้อหลัง MT

• 1. Function Approximation (1 ข ้อ)

– Taylor Series/McLauren Series

• 2. Roots of Function (1 ข ้อ)

– Bisection

– Newton-Ralphson

• 3. Linear Equations (1 ข ้อ) เรื่องใดเรื่องหนึ่ง

– Gauss Elimination

– Gauss Jordan (Including Matrix Inverse)

– Gauss Seidel

– LU Decomposition (Crout Decomposition)

Course Ends

• Prepare For Exam

• 4. Numerical Integration (1 ข ้อ) ไม่ออก Finite

Difference

– Trapezoidal Rule

– Simpson 1/3 Rule

– Richardson Extrapolation

• 5. ODE (1 ข ้อ)

– Classical Forth Order RK Method เรื่องเดียว

• 6. Curve Fitting (1 ข ้อ)

– Natural Spline

– OLS Regression

Formulas

ดูใน MT+ต่อไปนี้

END OF CPE 332 T1-57

Minimum 40% เพื่อที่จะผ่านวิชานี้

A จะต้องได้ 80% ขึ้นไป

Document Outline

Table of contents

previous page start